Method for determining a performance of an unsaturated subgrade of a multi-layered expressway section and processing device thereof

ABSTRACT

A method for determining a performance of a subgrade of an expressway section and a processing device thereof are used to accurately calculate a stress, a displacement, a velocity and an acceleration at any position inside the subgrade while analyzing the performance of the subgrade of the target expressway section, in full consideration of the moving speed and vibration characteristics of the vehicle load, the layered characteristics of the expressway and the unsaturated characteristics of the subgrade soil.

CROSS-REFERENCE TO RELATED APPLICATION

The present disclosure claims priority to Chinese Patent Application No. 202111460692.4 filed on Dec. 2, 2021, the disclosure of which is hereby incorporated by reference in its entirety.

BACKGROUND

An expressway is generally considered to be composed of a surface layer, a base layer and a subgrade. Under the action of vehicle load, the subgrade will bear a vibration load transmitted from an upper structure. Long-term vibration load during the service life will render the plastic deformation of subgrade soil gradually accumulated to thereby produce uneven settlement, which may lead to damage, like cavitation or cracking, of the upper surface structure. Thus, in order to guarantee long service life of expressways, it is necessary to quantify the mechanical properties of the subgrade in the road structure design.

At present, a commonly used calculation method is to simplify the structure of the expressway as a layered elastic medium, and simplify the vehicle load as a static load that is constant or vibrates at a certain frequency.

However, during the researches on current related technologies, the inventors found that in the real service state, the vehicle load is movable, and the subgrade soil is usually in an unsaturated state. Three phases, i.e., solid phase, liquid phase and gas phase, coexist in soil, and water and gas in the pores will greatly affect the mechanical properties of the soil. The layered elastic model completely neglects the multiple phases of soil, and does not take into consideration the influence of the change in water-containing state inside the subgrade on the entire structure, which will result in significant errors in the predication of the mechanical properties of the subgrade.

SUMMARY

The embodiments of the disclosure relate to the field of geology, and particularly to a method for determining a performance of an unsaturated subgrade of a multi-layered expressway section and a processing device thereof.

The present disclosure provides a method for determining a performance of an unsaturated subgrade of a multi-layered expressway section, as well as a processing device thereof, which are used to accurately calculate a stress, a displacement, a velocity and an acceleration at any position inside the subgrade while analyzing the performance of the subgrade of the target expressway section, in full consideration of the moving speed and vibration characteristics of the vehicle load, the layered characteristics of the expressway and the unsaturated characteristics of the subgrade soil.

In a first aspect, the present disclosure provides a method for determining a performance of an unsaturated subgrade of a multi-layered expressway section. The method may include:

Step S101: a processing device using an Odemark equivalent hypothesis of modulus-thickness to simplify a road structure of a target expressway section into a three-layered model, wherein the three-layered model is

${H_{e} = {\left( \frac{E}{E_{1}} \right)^{1/3}H}},$ wherein H_(e) is an equivalent thickness, E is a modulus of a layer to be transformed, H is a thickness of the layer to be transformed, and E₁ is a modulus of a target layer;

Step S102: based on the three-layered model, the processing device setting a base layer of the target expressway section as an elastic medium and setting the subgrade of the target expressway section as an unsaturated porous medium to derive and obtain a rigid pavement governing equation, a flexible pavement and base layer governing equation, and an unsaturated subgrade governing equation, wherein the rigid pavement governing equation is:

${{{D_{p}\left\lbrack {\frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial x^{4}} + {2\frac{\partial^{4}{W\left( {x,y,t} \right)}}{{\partial x^{2}}{\partial y^{2}}}} + \frac{\partial^{4}{w\left( {x,y,t} \right)}}{\partial y^{4}}} \right\rbrack} + {\rho_{p}h_{p}{\overset{¨}{W}\left( {x,y,t} \right)}} + {\varphi\left( {x,y,t} \right)}} = {f\left( {x,y,t} \right)}},$

wherein D_(p)=E_(p)h_(p) ³/[12(1−μ_(p) ²)] is a bending rigidity of a pavement slab, W(x, y, t) is a vertical displacement of the pavement slab, E_(p), ρ_(p) and h_(p) are respectively a modulus, a density and a thickness of the pavement slab, f(x, y, t) is a traffic load on the pavement slab, and φ(x, y, t) is a reaction force exerted by the base layer on the pavement slab,

the flexible pavement and base layer governing equation is

$\left\{ {\begin{matrix} {{{G_{b}u_{i,{jj}}^{b}} + {\left( {\lambda_{b} + G_{b}} \right)u_{j,{ji}}^{b}}} = {\rho_{b}{\overset{¨}{u}}_{i}^{b}}} \\ {\tau_{ij}^{b} = {{\lambda_{b}\theta_{b}\delta_{ij}} = {2G_{b}\varepsilon_{ij}^{b}}}} \end{matrix},} \right.$

wherein u^(b) is a displacement of an elastic medium layer, G_(b) and λ_(b) are Lame constants of the elastic medium layer, ρ_(b) is a density of the elastic medium layer, τ^(b) is a stress of the elastic medium layer, τ^(b) and ε^(b) respectively represent the stress of the elastic medium layer and a strain of the elastic medium layer, and δ_(ij) is a Kronecker symbol,

the unsaturated subgrade governing equation is

$\left\{ \begin{matrix} {{- p_{w,i}} = {{\rho_{w}{\overset{¨}{u}}_{i}} + {\frac{\rho_{w}}{{nS}_{r}}{\overset{¨}{W}}_{i}} + {\frac{\rho_{w}g}{k_{w}}{\overset{.}{W}}_{i}}}} \\ {{- p_{a,i}} = {{\rho_{a}{\overset{¨}{u}}_{i}} + {\frac{\rho_{a}}{n\left( {1 - S_{r}} \right)}{\overset{¨}{V}}_{i}} + {\frac{\rho_{a}g}{k_{a}}{\overset{.}{V}}_{i}}}} \\ {{{A_{11}{\overset{.}{p}}_{w}} + {A_{12}{\overset{.}{p}}_{a}} + {A_{13}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{14}{\nabla \cdot {\overset{.}{V}}_{i}}}} = 0} \\ {{{A_{21}{\overset{.}{p}}_{w}} + {A_{22}{\overset{.}{p}}_{a}} + {A_{23}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{24}{\nabla \cdot {\overset{.}{W}}_{i}}}} = 0} \\ {{{\mu{\nabla^{2}{\overset{.}{u}}_{i}}} + {\left( {\lambda + G} \right){\nabla e}} - {a\chi{\nabla p_{w}}} - {{a\left( {1 - \chi} \right)}{\nabla p_{a}}}} = {{\rho{\overset{¨}{u}}_{i}} + {\rho_{w}{\overset{¨}{W}}_{i}} + {\rho_{a}{\overset{¨}{V}}_{i}}}} \end{matrix} \right.$

wherein

${A_{11} = {{nA}_{s} + {\left( {1 - S_{r}} \right)\frac{{a\chi} - {nS}_{r}}{K_{s}}}}},$ ${A_{12} = {{\left( {1 - S_{r}} \right)\left\lbrack {\frac{{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}}{K_{s}} + \frac{n}{K_{a}}} \right\rbrack} - {nA}_{s}}},$ ${A_{13} = {\left( {1 - S_{r}} \right)\left( {1 - \frac{K_{b}}{K_{s}}} \right)}},{A_{14} = 1},$ ${A_{21} = {\frac{\left( {{a\chi} - {nS}_{r}} \right)S_{r}}{K_{s}} - {nA}_{s} + \frac{{nS}_{r}}{K_{w}}}},$ ${A_{22} = {\frac{\left\lbrack {{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}} \right\rbrack S_{r}}{K_{s}} + {nA}_{s}}},$ ${A_{23} = {\left( {1 - \frac{K_{b}}{K_{s}}} \right)S_{r}}},{A_{24} = 1},$ ${A_{s} = {{- \alpha}{{md}\left( {1 - S_{0}} \right)}{\left( S_{e} \right)^{\frac{m + 1}{m}}\left\lbrack {\left( S_{e} \right)^{- \frac{1}{m}} - 1} \right\rbrack}^{\frac{d - 1}{d}}}},$ a = 1 − K_(b)/K_(s), K_(b)andK_(s) are bulk compressibility moduli of soil skeletons and soil particles respectively, p_(w) and p_(a) are pressures of pore water and gas respectively, ρ_(w) and ρ_(a) are densities of pore water and gas respectively, W_(i) and V_(i) are displacement components of pore water and pore gas relative to soil particles in an i direction, g represents a gravitational acceleration, k_(w) and k_(a) represent permeability coefficients of pore water and pore gas respectively, χ is an effective stress parameter, α, m and d are fitting parameters of a soil-water characteristic curve model, S_(e) is an effective saturation, S_(w0) is an irreducible saturation, n is a soil porosity, and S_(r) is a saturation;

Step S103: based on an action form of a vehicle load, the processing device respectively setting time terms of variables in the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as a simple harmonic form e^(iΩ) ⁰ ^(t) and separating the simple harmonic form e^(iΩ) ⁰ ^(t) of the time terms from the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation, and then using a double Fourier transform to respectively reduce partial derivatives in a three-dimensional space domain of the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as ordinary differential equations, to obtain element stiffness matrices of different layer positions after solving and then establish a total dynamic stiffness matrix, wherein the double Fourier transform is

${{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)} = {\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{f\left( {x,y,z,t} \right)}e^{- {i({{\beta x} + {\gamma y}})}}{dxdy}}}}},$ ${{f\left( {x,y,z,t} \right)} = {\frac{1}{4\pi^{2}}{\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)}e^{i({{\beta x} + {\gamma y}})}d\beta d\gamma}}}}},$

wherein β and γ are wave numbers in x and y directions respectively, and f and f are variables in the space domain and a corresponding transform domain respectively,

the rigid pavement governing equation in the transform domain is: W (β,γ)[D _(p)(β²+γ²)²−Ω²ρ_(p) h _(p)]+φ(β,γ)= f (β,γ),

a stress-strain relationship of the flexible pavement in the transform domain is:

${\left\{ {\overset{¯}{\sum}}_{P} \right\} = {{{\left\lbrack S^{p} \right\rbrack\left\lbrack D^{p} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{P} \right\}} = {\left\lbrack K^{p} \right\rbrack\left\{ {\overset{¯}{U}}_{P} \right\}}}},{\left\lbrack D^{P} \right\rbrack = \begin{bmatrix} \frac{{- i}\beta{e^{- a_{1}}}^{h}}{k_{1}^{2}} & \frac{{- i}\beta}{k_{1}^{2}} & e^{{- \alpha_{2}}h} & 1 & 0 & 0 \\ \frac{{- i}\gamma{e^{- a_{1}}}^{h}}{k_{1}^{2}} & \frac{{- i}\gamma}{k_{1}^{2}} & 0 & 0 & e^{{- \alpha_{2}}h} & 1 \\ \frac{{- i}\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{i\alpha_{1}}{k_{1}^{2}} & {\frac{\beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \beta}{\alpha_{2}} & {\frac{\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \gamma}{\alpha_{2}} \\ \frac{{- i}\beta}{k_{1}^{2}} & \frac{{- i}\beta e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 1 & e^{{- \alpha_{2}}h} & 0 & 0 \\ \frac{{- i}\gamma}{k_{1}^{2}} & \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 0 & 0 & 1 & e^{{- \alpha_{2}}h} \\ \frac{{- i}\alpha_{1}}{k_{1}^{2}} & \frac{i\alpha_{1}{e^{- a_{1}}}^{h}}{k_{1}^{2}} & \frac{\beta}{\alpha_{2}} & {\frac{- \beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\gamma}{\alpha_{2}} & {\frac{- \gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \end{bmatrix}},{\left\lbrack S^{P} \right\rbrack = {G_{p}\begin{bmatrix} {\frac{2i\beta\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{2i\beta\alpha_{1}}{k_{1}^{2}} & {\frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{\alpha_{2}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{\beta\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} \\ {\frac{2i{\gamma\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{2i{\gamma\alpha}_{1}}{k_{1}^{2}} & {\frac{\beta\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {\frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} \\ {\frac{i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}} & {{- 2}\beta e^{{- \alpha_{2}}h}} & {{- 2}\beta} & {{- 2}\gamma e^{{- \alpha_{2}}h}} & {{- 2}\gamma} \\ \frac{{- 2}i\beta\alpha_{1}}{k_{1}^{2}} & {\frac{2i\beta\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {\frac{\beta\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{{- 2}i{\gamma\alpha}_{1}}{k_{1}^{2}} & {\frac{2i{\gamma\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {\frac{\beta\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}} & {\frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {2\beta} & {2\beta e^{{- \alpha_{2}}h}} & {2\gamma} & {2\gamma e^{{- \alpha_{2}}h}} \end{bmatrix}}},$

wherein superscript 0 of elements in the matrix represents a top surface of a flexible pavement layer, and superscript 1 of elements in the matrix represents a bottom surface of the flexible pavement layer, a stress vector

${\left\{ {\sum\limits^{\_}}_{P} \right\} = \left\{ {{- {\overset{¯}{\tau}}_{{xz},p}^{0}} - {\overset{¯}{\tau}}_{{yz},p}^{0} - {i{\overset{¯}{\tau}}_{{zz},p}^{0}{\overset{¯}{\tau}}_{{xz},p}^{1}{\overset{¯}{\tau}}_{{yz},p}^{1}i{\overset{¯}{\tau}}_{{zz},p}^{1}}} \right\}^{T}},\left\lbrack D^{p} \right\rbrack$ is a displacement matrix of the flexible pavement, [S^(p)] is a stress matrix of the flexible pavement, a displacement vector {Ū_(P)}={ū_(x,p) ⁰ ū_(y,p) ⁰ iū_(z,p) ⁰ ū_(x,p) ¹ ū_(y,p) ¹ iū_(z,p) ¹}^(T), [K^(p)]=[S^(p)][D^(p)]⁻¹ is a dynamic stiffness matrix of the flexible pavement layer, α₁ ²=β²+γ²−k₁ ², α₂ ²=β²+γ²−k₂ ², k₁ ²=ω²/c₁ ², k₂ ²=ω²/c₂ ², c is a traffic load velocity, c₁ and c₂ are a compression wave velocity and a shear wave velocity of the elastic medium, respectively, ω=Ω₀−βc, h represents a thickness of the flexible pavement layer, and i is an imaginary number,

a stress-strain relationship of the elastic base layer in the transform domain is: {Σ _(b) }=[S ^(b) ][D ^(b)]⁻¹ {Ū _(b) }=[K ^(b) ]{Ū _(b)},

a stress-strain relationship of the unsaturated subgrade in the transform domain is:

${\left\{ {\sum\limits^{\_}}_{s} \right\} = {{{\left\lbrack S^{s} \right\rbrack\left\lbrack D^{s} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{s} \right\}} = {\left\lbrack K^{s} \right\rbrack\left\{ {\overset{¯}{U}}_{s} \right\}}}},{\left\lbrack D^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma}{\beta}}e^{{- \lambda_{0}}z}} & {\frac{i\lambda_{0}}{\beta}e^{{- \lambda_{0}}z}} & {S_{1}i\beta e^{{- \lambda_{1}}z}} & {S_{2}i\beta e^{{- \lambda_{2}}z}} & {S_{3}i\beta e^{{- \lambda_{3}}z}} \\ e^{{- \lambda_{0}}z} & 0 & {S_{1}i\gamma e^{{- \lambda_{1}}z}} & {S_{2}i\gamma e^{{- \lambda_{2}}z}} & {S_{3}i\gamma e^{{- \lambda_{3}}z}} \\ 0 & {{- i}e^{{- \lambda_{0}}z}} & {{- i}S_{1}\lambda_{1}e^{{- \lambda_{1}}z}} & {{- i}S_{2}\lambda_{2}e^{{- \lambda_{2}}z}} & {{- i}S_{3}\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{w}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{w1}}{b^{w}\rho_{w}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{w2}}{b^{w}\rho_{w}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{w3}}{b^{w}\rho_{w}}} - e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{a}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{a1}}{b^{a}\rho_{a}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{a2}}{b^{a}\rho_{a}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{a3}}{b^{a}\rho_{a}}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$

${\left\lbrack S^{s} \right\rbrack = \begin{bmatrix} {{- \frac{{\gamma\lambda}_{0}G}{\beta}}e^{{- \lambda_{0}}z}} & {i{G\left( {\frac{\lambda_{0}^{2}}{\beta} + \beta} \right)}e^{{- \lambda_{0}}z}} & {2GS_{1}i\beta\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i\beta\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i\beta\lambda_{3}e^{{- \lambda_{3}}z}} \\ {G\lambda_{0}e^{{- \lambda_{0}}z}} & {{iG}\gamma e^{{- \lambda_{0}}z}} & {2GS_{1}i\gamma\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i{\gamma\lambda}_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i{\gamma\lambda}_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- i}2G\lambda_{0}e^{{- \lambda_{0}}z}} & {{- i}B_{1}e^{{- \lambda_{1}}z}} & {{- i}B_{2}e^{{- \lambda_{2}}z}} & {{- i}B_{3}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{w1}}e^{{- \lambda_{1}}z}} & {{- f_{w2}}e^{{- \lambda_{2}}z}} & {{- f_{w3}}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{a1}}e^{{- \lambda_{1}}z}} & {{- f_{a2}}e^{{- \lambda_{2}}z}} & {{- f_{a3}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$

wherein

$\begin{matrix} {{\left\{ \sum_{s} \right\} = \left\{ {{- {\overset{¯}{\tau}}_{{xz},s}}\  - {\overset{¯}{\tau}}_{{yz},s}\  - {i{\overset{¯}{\tau}}_{{zz},s}}\  - {\overset{¯}{p}}_{w}\  - {\overset{¯}{p}}_{a}} \right\}^{T}},} \\ \begin{matrix} {{\left\{ {\overset{¯}{U}}_{s} \right\} = \left\{ {{\overset{¯}{u}}_{x,s}\ {\overset{¯}{u}}_{y,s}\ i{\overset{¯}{u}}_{z,s}\ \overset{¯}{W}\ \overset{¯}{V}} \right\}^{T}},} & {{B_{wn} = {{\omega^{2}\rho_{w}S_{n}\lambda_{n}} - {f_{wn}\lambda_{n}}}},} \end{matrix} \\ \begin{matrix} {{B_{an} = {{\omega^{2}\rho_{a}S_{n}\lambda_{n}} - {f_{an}\lambda_{n}}}},} & {{B_{n} = {{2\mu S_{n}\lambda_{n}^{2}} + \lambda - {a\chi f_{wn}} - {a\left( {1 - \chi} \right)f_{an}}}},} \end{matrix} \\ \begin{matrix} {{S_{n} = {- \frac{\lambda + \mu + {b_{1}f_{wn}} + {b_{2}f_{an}}}{b_{3} + {\mu d_{n}}}}},} & {{f_{an} = \frac{{b_{13}b_{31}} - {b_{11}b_{33}} + {b_{11}d_{n}}}{{b_{11}b_{32}} - {b_{31}b_{12}} + {b_{31}d_{n}}}},} & {{f_{m} = \frac{{b_{32}b_{23}} - {b_{22}b_{33}} + {b_{22}d_{n}}}{{b_{22}b_{31}} - {b_{32}b_{21}} + {b_{32}d_{n}}}},} \end{matrix} \\ {\begin{matrix} {{\lambda_{n} = \sqrt{\beta^{2} + \gamma^{2} + d_{n}}},} & {{d_{n}{is}{roots}{of}{an}{equation}^{{d_{n}^{3} + {b_{4}d_{n}^{2}} + {b_{5}d_{n}} + b_{6}} = 0}},} & {b_{4} = {{- b_{12}} - b_{33} - b_{21}}} \end{matrix},} \\ {b_{5} = {{b_{12}b_{33}} - {b_{13}b_{32}} + {b_{21}b_{12}} + {b_{21}b_{33}} - {b_{22}b_{11}} - {b_{23}b_{31}}}} \\ \begin{matrix} {{b_{6} = {{{- b_{21}}b_{12}b_{33}} - {b_{13}b_{31}b_{22}} - {b_{23}b_{11}b_{32}} + {b_{21}b_{13}b_{32}} + {b_{22}b_{11}b_{33}} + {b_{23}b_{31}b_{12}}}},} & {{n = 1},2,3,} & {{b_{1} = {{{- a}\chi} - \frac{\omega^{2}}{b^{w}}}},} \end{matrix} \\ \begin{matrix} {{b_{2} = {{- {a\left( {1 - \chi} \right)}} - \frac{\omega^{2}}{b^{a}}}},} & {{b_{3} = {{\rho\omega}^{2} + \frac{\omega^{4}\rho_{w}}{b^{w}} + \frac{\omega^{4}\rho_{a}}{b^{a}}}},} & {{b_{11} = {b^{a}\rho_{a}A_{11}}},} & {b_{12} = {b^{a}\rho_{a}A_{12}}} \end{matrix} \\ \begin{matrix} {{b_{13} = {{\omega^{2}\rho_{a}} + {b^{a}\rho_{a}A_{13}}}},} & {{b_{21} = {b^{w}\rho_{w}A_{21}}},} & {{b_{22} = {b^{w}\rho_{w}A_{22}}},} & {{b_{23} = {{\omega^{2}\rho_{w}} + {b^{w}\rho_{w}A_{23}}}},} \end{matrix} \\ {\begin{matrix} {{b_{31} = {- \frac{{b_{1}b_{21}} + {b_{2}b_{11}}}{\lambda + {2\mu}}}},} & {{b_{32} = {- \frac{{b_{1}b_{22}} + {b_{2}b_{12}}}{\lambda + {2\mu}}}},} & {{b_{33} = {- \frac{{b_{1}b_{23}} + {b_{2}b_{13}} + b_{3}}{\lambda + {2\mu}}}},} & {b^{w} = {{- \frac{\omega^{2}}{nS_{r}}} + \frac{i\omega g}{k_{w}}}} \end{matrix},} \\ {{b^{a} = {{- \frac{\omega^{2}}{n\left( {1 - S_{r}} \right)}} + \frac{i\omega g}{k_{a}}}},} \end{matrix}$

a total stiffness matrix of a flexible pavement expressway model is: {Σ _(f) }=[K ^(f) ]{Ū _(f)},

wherein {Σ _(f)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ 0 0 0 0 0 0 0 0}^(T), {Ū_(f)}={ū_(x) ⁰ ū_(y) ⁰ iū_(z) ⁰ ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), superscripts 0, 1, 2 correspond to different layer positions shown from top to bottom in the three-layered model, and [K^(f)] is a matrix of 11×11,

a total stiffness matrix of a rigid pavement expressway model is: {Σ _(r) }=[K ^(r) {Ū _(r)},

wherein {Σ _(r)}={−τ _(xz,b) ¹ −τ _(yz,b) ¹ −iτ _(zz,b) ¹ 0 0 0 0 0}^(T), {Ū_(r)}={ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), and [K^(r)] is a matrix of 8×8;

Step S104: the processing device determining a load form and a load size according to a specific vehicle, substituting the load form and the load size as boundary conditions into the stiffness matrix of a corresponding model to solve and obtain a solution f of a parameter to be solved in a frequency domain, wherein the load form and the load size determined based on the specific vehicle are

${f\left( {x_{1},y_{1},t} \right)} = \left\{ \begin{matrix} {\frac{Pe^{i\Omega_{0}t}}{40{l}_{1}l_{2}};\left\{ \begin{matrix} {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + l_{2}} \right)}} \end{matrix} \right.} \\ {0;{others}} \end{matrix} \right.$

wherein P represents a total axle load of the specific vehicle, and in a moving coordinate system, x=x₁−Vt, y=y₁, z=z₁;

wherein in the transform domain, the load form and the load size determined based on the specific vehicle are

${{\overset{¯}{f}\left( {\beta,\gamma} \right)} = \frac{P\left( {e^{{- i}\gamma d_{w}} + e^{y\gamma d_{w}}} \right)}{10l_{1}l_{2}\beta\gamma}}\text{ }\left\lbrack {{{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}e^{{- i}\beta d_{c}}} + {{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}\left( {e^{{- i}\beta d_{a}} + e^{i\beta d_{a}}} \right)\left( {e^{{- i}\gamma d_{t}}\  + e^{i\gamma d_{t}}} \right)}} \right\rbrack$

a boundary condition of the rigid pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},b}^{1}\left( {x,y,t} \right)} = {- {\varphi\left( {x,\ y,t} \right)}}} \\ {{\tau_{{xz},b}^{1}\left( {x,y,t} \right)} = {{\tau_{{yz},b}^{1}{p\left( {x,y,t} \right)}} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \\ {{u_{{zz},b}^{1}\left( {x,y,t} \right)} = {W\left( {x,y,t} \right)}} \end{matrix},} \right.$

the boundary condition of the rigid pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {- {\overset{¯}{\varphi}\left( {\beta,\gamma} \right)}}} \\ {{{\overset{¯}{\tau}}_{{xz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{{yz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \\ {{{\overset{¯}{u}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {\overset{\_}{W}\left( {\beta,\gamma} \right)}} \end{matrix},} \right.$

wherein

${{\overset{¯}{\varphi}\left( {\beta,\gamma} \right)} = \frac{\overset{¯}{f}\left( {\beta,\gamma} \right)}{1 + {B_{33}{D_{P}\left( {\beta^{2} + \gamma^{2}} \right)}^{2}} - {B_{33}\omega^{2}m_{b}}}},$ and B₃₃ is an element in a third row and a third column of [K^(r)]⁻¹ and corresponds to a vertical displacement of the base layer,

a boundary condition of the flexible pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},p}^{0}\left( {x,y,t} \right)} = {f\left( {x,y,t} \right)}} \\ {{\tau_{xz}^{0}\left( {x,y,t} \right)} = {{\tau_{yz}^{0}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \end{matrix},} \right.$

the boundary condition of the flexible pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},p}^{0}\left( {\beta,\gamma} \right)} = {\overset{¯}{f}\left( {\beta,\gamma} \right)}} \\ {{{\overset{¯}{\tau}}_{xz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{yz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \end{matrix},} \right.$

wherein for the rigid pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the rigid pavement expressway model, and for the flexible pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the flexible pavement expressway model, in such a way to obtain the displacement of each surface, the displacement {ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²} of the top surface of the subgrade is substituted into the stress-strain relationship of the unsaturated subgrade in the transform domain to obtain the a stress {τ _(xz,s) τ _(yz,s) τ _(zz,s) p _(w) p _(α)} and a displacement {ū_(x) ² ū_(y) ² ū_(z) ² W ² V ²} at any position inside the subgrade in the frequency domain, and a vibration velocity v _(ij)=iωū_(ij) and an acceleration ā_(ij)=−ωū_(ij) of the subgrade can be calculated according to the displacement;

Step S105: the processing device using an inverse fast Fourier transform to convert the solution f of the parameter to be solved in the frequency domain into a solution f in the time domain, wherein the solution f in the time domain is a stress {τ_(xz,s) τ_(yz,s) τ_(zz,s) p_(w) p_(α)}, a displacement {u² _(x) u² _(y) u² _(z) W² V²}, a velocity v_(ij) and an acceleration a_(ij).

In combination of the first aspect of the present disclosure, in the first possible implementation mode according to the first aspect of the present disclosure, the method further includes Step S106, which is prior to Step S101:

the processing device acquiring a task of determining the performance of the subgrade of the target expressway section.

In combination of the first possible implementation mode according to the first aspect of the present disclosure, in the second possible implementation mode according to the first aspect of the present disclosure, Step S106 specifically includes:

the processing device monitoring whether project data on a system are updated; and

in response to monitoring that the project data are updated, the processing device acquiring the task of determining the performance of the subgrade of the target expressway section according to the updated project data.

In combination of the first possible implementation mode according to the first aspect of the present disclosure, in the third possible implementation mode according to the first aspect of the present disclosure, Step S106 specifically includes:

the processing device receiving the task of determining the performance of the subgrade of the target expressway section as inputted by a user.

In combination of the first aspect of the present disclosure, in the fourth possible implementation mode according to the first aspect of the present disclosure, the method further includes Step S107, which is subsequent to Step S105:

the processing device outputting the solution f in the time domain.

In combination of the fourth possible implementation mode according to the first aspect of the present disclosure, in the fifth possible implementation mode according to the first aspect of the present disclosure, Step S107 specifically includes:

the processing device outputting a performance determination report of the subgrade of the target expressway section, wherein the performance determination report includes the solution f in the time domain.

In combination of the first aspect of the present disclosure, in the sixth possible implementation mode according to the first aspect of the present disclosure, the method further includes Step S108, which is subsequent to Step S105:

the processing device determining whether a value of the solution f in the time domain reaches a predetermined alarm value; and

in response to determining the value of the solution f in the time domain reaches the predetermined alarm value, the processing device outputting an alarm prompt.

In a second aspect, the present disclosure provides a processing device for determining a performance of an unsaturated subgrade of a multi-layered expressway section. The processing device includes a processor and a memory in which a computer program is stored, wherein when the computer program in the memory is executed, the processor is configured to:

use an Odemark equivalent hypothesis of modulus-thickness to simplify a road structure of a target expressway section into a three-layered model, wherein the three-layered model is

${H_{e} = {\left( \frac{E}{E_{1}} \right)^{1/3}H}},$ wherein H_(e) is an equivalent thickness, E is a modulus of a layer to be transformed, H is a thickness of the layer to be transformed, and E₁ is a modulus of a target layer;

based on the three-layered model, set a base layer of the target expressway section as an elastic medium and set the subgrade of the target expressway section as an unsaturated porous medium to derive and obtain a rigid pavement governing equation, a flexible pavement and base layer governing equation, and an unsaturated subgrade governing equation, wherein the rigid pavement governing equation is:

${{{D_{p}\left\lbrack {\frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial x^{4}} + {2\frac{\partial^{4}{W\left( {x,y,t} \right)}}{{\partial x^{3}}{\partial y^{2}}}} + \frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial x^{4}}} \right\rbrack} + {\rho_{p}h_{p}{\overset{¨}{W}\left( {x,y,t} \right)}} + {\varphi\left( {x,y,t} \right)}} = {f\left( {x,y,t} \right)}},$

wherein D_(p)=E_(p)h_(p) ³/[12(1−μ_(p) ²)] is a bending rigidity of a pavement slab, W(x, y, t) is a vertical displacement of the pavement slab, E_(p), ρ_(p) and h_(p) are respectively a modulus, a density and a thickness of the pavement slab, f(x, y, t) is a traffic load on the pavement slab, and φ(x, y, t) is a reaction force exerted by the base layer on the pavement slab,

the flexible pavement and base layer governing equation is

$\left\{ {\begin{matrix} {{{G_{b}u_{i,{jj}}^{b}} + {\left( {\lambda_{b} + G_{b}} \right)u_{j,{ji}}^{b}}} = {\rho_{b}{\overset{¨}{u}}_{i}^{b}}} \\ {\tau_{ij}^{b} = {{\lambda_{b}\theta_{b}\delta_{ij}} + {2G_{b}\varepsilon_{ij}^{b}}}} \end{matrix},} \right.$

wherein u^(b) is a displacement of an elastic medium layer, G_(b) and λ_(b) are Lame constants of the elastic medium layer, ρ_(b) is a density of the elastic medium layer, τ^(b) is a stress of the elastic medium layer, τ^(b) and ε^(b) respectively represent the stress of the elastic medium layer and a strain of the elastic medium layer, and δ_(ij) is a Kronecker symbol,

the unsaturated subgrade governing equation is

$\left\{ {\begin{matrix} {{- p_{w,i}} = {{\rho_{w}{\overset{¨}{u}}_{i}} + {\frac{\rho_{w}}{{nS}_{r}}{\overset{¨}{W}}_{i}} + {\frac{\rho_{w}g}{k_{w}}{\overset{.}{W}}_{i}}}} \\ {{- p_{a,i}} = {{\rho_{a}{\overset{¨}{u}}_{i}} + {\frac{\rho_{a}}{n\left( {1 - S_{r}} \right)}{\overset{¨}{V}}_{i}} + {\frac{\rho_{a}g}{k_{a}}{\overset{.}{V}}_{i}}}} \\ {{{A_{11}{\overset{.}{p}}_{w}} + {A_{12}{\overset{.}{p}}_{a}} + {A_{13}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{14}{\nabla \cdot {\overset{.}{V}}_{i}}}} = 0} \\ {{{A_{21}{\overset{.}{p}}_{w}} + {A_{22}{\overset{.}{p}}_{a}} + {A_{23}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{24}{\nabla \cdot {\overset{.}{W}}_{i}}}} = 0} \\ {{{\mu{\nabla^{2}{\overset{.}{u}}_{i}}} + {\left( {\lambda + G} \right){\nabla e}} - {a\chi{\nabla p_{w}}} - {{a\left( {1 - \chi} \right)}{\nabla p_{a}}}} = {{\rho{\overset{¨}{u}}_{i}} + {\rho_{w}{\overset{¨}{W}}_{t}} + {\rho_{a}{\overset{¨}{V}}_{i}}}} \end{matrix},} \right.$

wherein

${A_{11} = {{nA_{s}} + {\left( {1 - S_{r}} \right)\frac{{a\chi} - {nS_{r}}}{K_{s}}}}},$ ${A_{12} = {{\left( {1 - S_{r}} \right)\left\lbrack {\frac{{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}}{K_{s}} + \frac{n}{K_{a}}} \right\rbrack} - {nA_{s}}}},$ ${A_{13} = {\left( {1 - S_{r}} \right)\left( {1 - \frac{K_{b}}{K_{s}}} \right)}},$ A₁₄ = 1, ${A_{21} = {\frac{\left( {{a\chi} - {nS_{r}}} \right)S_{r}}{K_{s}} - {nA_{s}} + \frac{nS_{r}}{K_{w}}}},$ ${A_{22} = {\frac{\left\lbrack {{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}} \right\rbrack S_{r}}{K_{s}} + {nA_{s}}}},$ ${A_{23} = {\left( {1 - \frac{K_{b}}{K_{s}}} \right)S_{r}}},$ A₂₄ = 1, ${A_{s} = {{- \alpha}{{md}\left( {1 - S_{w0}} \right)}{\left( S_{e} \right)^{\frac{m + 1}{m}}\left\lbrack {\left( S_{e} \right)^{- \frac{1}{m}} - 1} \right\rbrack}^{\frac{d - 1}{d}}}},$ a = 1 − K_(b)/K_(s), compressibility moduli of soil skeletons and soil particles respectively, p_(w) and p_(a) are pressures of pore water and gas respectively, ρ_(w) and ρ_(a) are densities of pore water and gas respectively, W_(i) and V_(i) are displacement components of pore water and pore gas relative to soil particles in an i direction, g represents a gravitational acceleration, k_(w) and k_(a) represent permeability coefficients of pore water and pore gas respectively, χ is an effective stress parameter, α, m and d are fitting parameters of a soil-water characteristic curve model, S_(e) is an effective saturation, S_(w0) is an irreducible saturation, n is a soil porosity, and S_(r) is a saturation;

wherein the processor is further configured to, based on an action form of a vehicle load, respectively set time terms of variables in the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation respectively as a simple harmonic form e^(iΩ) ⁰ ^(t) and separate the simple harmonic form e^(iΩ) ⁰ ^(t) of the time terms from the rigid pavement governing equation, the flexible pavement and base layer governing equation, and the unsaturated subgrade governing equation, and then use a double Fourier transform to respectively reduce partial derivatives in a three-dimensional space domain of the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as ordinary differential equations, to obtain element stiffness matrices of different layer positions after solving and then establish a total dynamic stiffness matrix, wherein the double Fourier transform is

${{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)} = {\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{f\left( {x,y,z,t} \right)}e^{- {({{\beta x} + {\gamma y}})}}{dxdy}}}}},$ ${{f\left( {x,y,z,t} \right)} = {\frac{1}{4\pi^{2}}{\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)}e^{i({{\beta x} + {\gamma y}})}d{\beta d\gamma}}}}}},$

wherein β and γ are wave numbers in x and y directions respectively, and f and f are variables in the space domain and a corresponding transform domain respectively,

the rigid pavement governing equation in the transform domain is: W (β,γ)[D _(p)(β²+γ²)²−Ω²ρ_(p) h _(p)]+φ(β,γ)= f (β,γ),

a stress-strain relationship of the flexible pavement in the transform domain is:

${\left\{ {\overset{¯}{\sum}}_{P} \right\} = {{{\left\lbrack S^{p} \right\rbrack\left\lbrack D^{p} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{P} \right\}} = {\left\lbrack K^{p} \right\rbrack\left\{ {\overset{¯}{U}}_{P} \right\}}}},$ [D^(P)]= $\begin{bmatrix} \frac{{- i}\beta e^{{- \alpha_{1}}h}}{{k_{1}}^{2}} & \frac{{- i}\beta}{{k_{1}}^{2}} & e^{{- \alpha_{2}}h} & 1 & 0 & 0 \\ \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{{k_{1}}^{2}} & \frac{{- i}\gamma}{{k_{1}}^{2}} & 0 & 0 & e^{{- \alpha_{2}}h} & 1 \\ \frac{{- {i\alpha}_{1}}e^{{- \alpha_{1}}h}}{{k_{1}}^{2}} & \frac{i\alpha_{1}}{{k_{1}}^{2}} & {\frac{\beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \beta}{\alpha_{2}} & {\frac{\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \gamma}{\alpha_{2}} \\ \frac{{- i}\beta}{{k_{1}}^{2}} & \frac{{- i}\beta e^{{- \alpha_{1}}h}}{{k_{1}}^{2}} & 1 & e^{{- \alpha_{2}}h} & 0 & 0 \\ \frac{{- i}\gamma}{{k_{1}}^{2}} & \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{{k_{1}}^{2}} & 0 & 0 & 1 & e^{{- \alpha_{2}}h} \\ \frac{{- i}\alpha_{1}}{{k_{1}}^{2}} & \frac{i\alpha_{1}e^{{- \alpha_{1}}h}}{{k_{1}}^{2}} & \frac{\beta}{\alpha_{2}} & {\frac{- \beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\gamma}{\alpha_{2}} & {\frac{- \gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \end{bmatrix},$

${\left\lbrack S^{P} \right\rbrack = {G_{p}\begin{bmatrix} {\frac{2i{\beta\alpha}_{1}}{{k_{1}}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i{\beta\alpha}_{1}}{{k_{1}}^{2}}} & {\frac{\left( {{\alpha_{2}}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\left( {{\alpha_{2}}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} \\ {\frac{2i{\gamma\alpha}_{1}}{{k_{1}}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i{\gamma\alpha}}{{k_{1}}^{2}}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\left( {{\alpha_{2}}^{2} + \gamma^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {{\alpha_{2}}^{2} + \gamma^{2}} \right)}{\alpha_{2}} \\ {\frac{i\left( {{2{\alpha_{2}}^{2}} + {k_{2}}^{2}} \right)}{{k_{1}}^{2}}e^{{- \alpha_{1}}h}} & \frac{i\left( {{2{\alpha_{2}}^{2}} + {k_{2}}^{2}} \right)}{{k_{1}}^{2}} & {{- 2}\beta e^{{- \alpha_{2}}h}} & {{- 2}\beta} & {{- 2}\gamma e^{{- \alpha_{2}}h}} & {{- 2}\gamma} \\ \frac{{- 2}i{\beta\alpha}_{1}}{{k_{1}}^{2}} & {\frac{2i{\beta\alpha}_{1}}{{k_{1}}^{2}}e^{{- \alpha_{1}}h}} & \frac{\left( {{\alpha_{2}}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {{\alpha_{2}}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} \\ \frac{{- 2}i\gamma\alpha_{1}}{{k_{1}}^{2}} & {\frac{2i\gamma\alpha}{{k_{1}}^{2}}e^{{- \alpha_{1}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {{\alpha_{2}}^{2} + \gamma^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {{\alpha_{2}}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{- {i\left( {{2{\alpha_{2}}^{2}} + {k_{2}}^{2}} \right)}}{{k_{1}}^{2}} & {\frac{- {i\left( {{2{\alpha_{2}}^{2}} + {k_{2}}^{2}} \right)}}{{k_{1}}^{2}}e^{{- \alpha_{1}}h}} & {2\beta} & {2\beta e^{{- \alpha_{2}}h}} & {2\gamma} & {2\gamma e^{{- \alpha_{2}}h}} \end{bmatrix}}},$

wherein superscript 0 of elements in the matrix represents a top surface of a flexible pavement layer, and superscript 1 of elements in the matrix represents a bottom surface of the flexible pavement layer, a stress vector {Σ _(P)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ τ _(xz,p) ¹ −τ _(yz,p) ¹ −iτ _(zz,p) ¹}^(T), [D^(p)] is a displacement matrix of the flexible pavement, [S^(p)] is a stress matrix of the flexible pavement, a displacement vector {Ū_(P)}={ū_(x,p) ⁰ ū_(y,p) ⁰ iū_(z,p) ⁰ ū_(x,p) ¹ ū_(y,p) ¹ iū_(z,p) ¹}^(T), [K^(p)]=[S^(p)][D^(p)]⁻¹ is a dynamic stiffness matrix of the flexible pavement layer, α₁ ²=β²+γ²−k₁ ², α₂ ²=β²+γ²−k₂ ², k₁ ²=ω²/c₁ ², k₂ ²=ω²/c₂ ², c is a traffic load velocity, c₁ and c₂ are a compression wave velocity and a shear wave velocity of the elastic medium respectively, ω=Ω₀−βc, h represents a thickness of the flexible pavement layer, and i is an imaginary number,

a stress-strain relationship of the elastic base layer in the transform domain is: {Σ _(b) }=[S ^(b) ][D ^(b)]⁻¹ {Ū _(b) }=[K ^(b) ]{Ū _(b)},

a stress-strain relationship of the unsaturated subgrade in the transform domain is:

${\left\{ {\sum\limits^{\_}}_{s} \right\} = {{{\left\lbrack S^{s} \right\rbrack\left\lbrack D^{s} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{s} \right\}} = {\left\lbrack K^{s} \right\rbrack\left\{ {\overset{¯}{U}}_{s} \right\}}}},$ ${\left\lbrack D^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma}{\beta}}e^{{- \lambda_{0}}z}} & {\frac{i\lambda_{0}}{\beta}e^{{- \lambda_{0}}z}} & {S_{1}i\beta e^{{- \lambda_{1}}z}} & {S_{2}i\beta e^{{- \lambda_{2}}z}} & {S_{3}i\beta e^{{- \lambda_{3}}z}} \\ e^{{- \lambda_{0}}z} & 0 & {S_{1}i\gamma e^{{- \lambda_{1}}z}} & {S_{2}i\gamma e^{{- \lambda_{2}}z}} & {S_{3}i\gamma e^{{- \lambda_{3}}z}} \\ 0 & {{- i}e^{{- \lambda_{0}}z}} & {{- i}S_{1}\lambda_{1}e^{{- \lambda_{1}}z}} & {{- i}S_{2}\lambda_{2}e^{{- \lambda_{2}}z}} & {{- i}S_{3}\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{w}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{w1}}{b^{w}\rho_{w}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{w2}}{b^{w}\rho_{w}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{w3}}{b^{w}\rho_{w}}}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{a}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{a1}}{b^{a}\rho_{a}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{a2}}{b^{a}\rho_{a}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{a3}}{b^{a}\rho_{a}}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$ ${\left\lbrack S^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma\lambda_{0}G}{\beta}}e^{{- \lambda_{0}}z}} & {i{G\left( {\frac{{\lambda_{0}}^{2}}{\beta} + \beta} \right)}e^{{- \lambda_{0}}z}} & {2GS_{1}i\beta\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i\beta\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i{\beta\lambda}_{3}e^{{- \lambda_{3}}z}} \\ {G\lambda_{0}e^{{- \lambda_{0}}z}} & {iG\gamma e^{{- \lambda_{0}}z}} & {2GS_{1}i\gamma\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i{\gamma\lambda}_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i{\gamma\lambda}_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- i}2G\lambda_{0}e^{{- \lambda_{0}}z}} & {{- i}B_{1}e^{{- \lambda_{1}}z}} & {{- i}B_{2}e^{{- \lambda_{2}}z}} & {{- i}B_{3}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{w1}}e^{{- \lambda_{1}}z}} & {{- f_{w2}}e^{{- \lambda_{2}}z}} & {{- f_{w3}}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{a1}}e^{{- \lambda_{1}}z}} & {{- f_{a2}}e^{{- \lambda_{2}}z}} & {{- f_{a3}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$

wherein

${\left\{ {\overset{-}{\sum}}_{s} \right\} = \left\{ {{- {\overset{¯}{\tau}}_{{xz},s}}\  - {\overset{¯}{\tau}}_{{yz},s}\  - {i{\overset{¯}{\tau}}_{{zz},s}}\  - {\overset{¯}{p}}_{w}\  - {\overset{¯}{p}}_{a}} \right\}^{T}},$ $\left\{ {{{\overset{\_}{U}}_{s} = \left\{ {{\overset{\_}{u}}_{x,s}{\overset{\_}{u}}_{y,s}i{\overset{\_}{u}}_{z,s}\overset{\_}{W}\overset{\_}{V}} \right\}^{T}},} \right.$ B_(wn) = ω²ρ_(w)S_(n)λ_(n) − f_(wn)λ_(n), B_(an) = ω²ρ_(a)S_(n)λ_(n) − f_(an)λ_(n), B_(n) = 2μS_(n)λ_(n)² + λ − aχf_(wn) − a(1 − χ)f_(an,) ${S_{n} = {- \frac{\lambda + \mu + {b_{1}f_{wn}} + {b_{2}f_{an}}}{b_{3} + {\mu d_{n}}}}},$ ${f_{an} = \frac{{b_{13}b_{31}} - {b_{11}b_{33}} + {b_{11}d_{n}}}{{b_{11}b_{32}} - {b_{31}b_{12}} + {b_{31}d_{n}}}},$ ${f_{wn} = \frac{{b_{32}b_{23}} - {b_{22}b_{33}} + {b_{22}d_{n}}}{{b_{22}b_{31}} - {b_{32}b_{21}} + {b_{32}d_{n}}}},$ ${\lambda_{n} = \sqrt{\beta^{2} + \gamma^{2} + d_{n}}},$ d_(n)isrootsofanequation d_(n)³ + b₄d_(n)² + b₅d_(n) + b₆ = 0, b₄ = −b₁₂ − b₃₃ − b₂₁, b₅ = b₁₂b₃₃ − b₁₃b₃₂ + b₂₁b₁₂ + b₂₁b₃₃ − b₂₂b₁₁ − b₂₃b₃₁, b₆ = −b₂₁b₁₂b₃₃ − b₁₃b₃₁b₂₂ − b₂₃b₁₁b₃₂ + b₂₁b₁₃b₃₂ + b₂₂b₁₁b₃₃ + b₂₃b₃₁b₁₂, n = 1, 2, 3, ${b_{1} = {{{- a}\chi} - \frac{\omega^{2}}{b^{w}}}},$ ${b_{2} = {{- {a\left( {1 - \chi} \right)}} - \frac{\omega^{2}}{b^{a}}}},$ ${b_{3} = {{\rho\omega}^{2} + \frac{\omega^{4}\rho_{w}}{b^{w}} + \frac{\omega^{4}\rho_{a}}{b^{a}}}},$ b₁₁ = b^(a)ρ_(a)A₁₁ b₁₂ = b^(a)ρ_(a)A₁₂, b₁₃ = ω²ρ_(a) + b^(a)ρ_(a)A₁₃, b₂₁ = b^(w)ρ_(w)A₂₁, b₂₂ = b^(w)ρ_(w)A₂₂, b₂₃ = ω²ρ_(w) + b^(w)ρ_(w)A₂₃, ${b_{31} = {- \frac{{b_{1}b_{21}} + {b_{2}b_{11}}}{\lambda + {2\mu}}}},$ ${b_{32} = {- \frac{{b_{1}b_{22}} + {b_{2}b_{12}}}{\lambda + {2\mu}}}},$ ${b_{33} = {- \frac{{b_{1}b_{23}} + {b_{2}b_{13}} + b_{3}}{\lambda + {2\mu}}}},$ ${b^{w} = {{- \frac{\omega^{2}}{{nS}_{r}}} + \frac{i\omega g}{k_{w}}}},$ ${b^{a} = {{- \frac{\omega^{2}}{n\left( {1 - S_{r}} \right)}} + \frac{i\omega g}{k_{a}}}},$

a total stiffness matrix of a flexible pavement expressway model is: {Σ _(f) }=[K _(f) ]{Ū _(f)},

wherein {Σ _(f)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ 0 0 0 0 0 0 0 0}^(T), {Ū_(f)}={ū_(x) ⁰ ū_(y) ⁰ iū_(z) ⁰ ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), superscripts 0, 1, 2 correspond to different layer positions shown from top to bottom in the three-layered model, and [K^(f)] is a matrix of 11×11,

a total stiffness matrix of a rigid pavement expressway model is: {Σ _(r) }=[K ^(r) ]{Ū _(r)},

wherein {Σ _(r)}={−τ _(xz,b) ¹ −τ _(yz,b) ¹ −iτ _(zz,b) ¹ 0 0 0 0 0}^(T), {Ū_(r)}={ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), and [K^(r)] is a matrix of 8×8;

wherein the processor is further configured to determine a load form and a load size according to a specific vehicle, substituting the load form and the load size as boundary conditions into the stiffness matrix of a corresponding model to solve and obtain a solution f of a parameter to be solved in a frequency domain, wherein the load form and the load size determined based on the specific vehicle are

${f\left( {x_{1},y_{1},t} \right)} = \left\{ {\begin{matrix} {\frac{{Pe}^{i\Omega_{a}x}}{40l_{1}l_{2}};\left\{ \begin{matrix} {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} + d_{1} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{1} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - d_{1} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{1} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} + d_{1} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{1} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} - d_{1} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{1} + l_{2}} \right)}} \\ {{\left( {d_{c} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{c} + l_{1}} \right)},{\left( {d_{w} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + l_{2}} \right)}} \end{matrix} \right.} \\ {0;{others}} \end{matrix},} \right.$

wherein P represents a total axle load of the specific vehicle, and in a moving coordinate system, x=x₁−Vt, y=y₁, z=z₁, wherein in the transform domain, the load form and the load size determined based on the specific vehicle are

${{\overset{¯}{f}\left( {\beta,\gamma} \right)} = {\frac{P\left( {e^{{- i}\gamma d_{w}} + e^{i\gamma d_{w}}} \right)}{10l_{1}l_{2}\beta\gamma}\left\lbrack {{{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}e^{{- i}\beta d_{c}}} + {{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}\left( {e^{{- i}\beta d_{a}} + e^{i\beta d_{a}}} \right)\left( {e^{{- i}\gamma d_{t}}\  + e^{i\gamma d_{t}}} \right.}} \right\rbrack}},$

a boundary condition of the rigid pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},b}^{1}\left( {x,y,t} \right)} = {- {\varphi\left( {x,y,t} \right)}}} \\ {{\tau_{{xz},b}^{1}\left( {x,y,t} \right)} = {{\tau_{{\gamma z},b}^{1}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \\ {{u_{{zz},b}^{1}\left( {x,y,t} \right)} = {W\left( {x,y,t} \right)}} \end{matrix},} \right.$

the boundary condition of the rigid pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{\_}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {- {\overset{\_}{\varphi}\left( {\beta,\gamma} \right)}}} \\ {{{\overset{\_}{\tau}}_{{xz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{\tau}}_{{\gamma z},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \\ {{{\overset{\_}{u}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {\overset{\_}{W}\left( {\beta,\gamma} \right)}} \end{matrix},} \right.$

wherein

${{\overset{¯}{\varphi}\left( {\beta,\gamma} \right)} = \frac{\overset{¯}{f}\left( {\beta,\gamma} \right)}{1 + {B_{33}{D_{P}\left( {\beta^{2} + \gamma^{2}} \right)}^{2}} - {B_{33}\omega^{2}m_{b}}}},$ and B₃₃ is an element in a third row and a third column of [K^(r)]⁻¹ and corresponds to a vertical displacement of the base layer,

a boundary condition of the flexible pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},p}^{0}\left( {x,y,t} \right)} = {f\left( {x,y,t} \right)}} \\ {{\tau_{xz}^{0}\left( {x,y,t} \right)} = {{\tau_{yz}^{0}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \end{matrix},} \right.$

the boundary condition of the flexible pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},p}^{0}\left( {\beta,\gamma} \right)} = {\overset{¯}{f}\left( {\beta,\gamma} \right)}} \\ {{{\overset{¯}{\tau}}_{xz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{yz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \end{matrix},} \right.$

wherein for the rigid pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the rigid pavement expressway model, and for the flexible pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the flexible pavement expressway model, in such a way to obtain the displacement of each surface, the displacement {ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²} of the top surface of the subgrade is substituted into the stress-strain relationship of the unsaturated subgrade in the transform domain to obtain a stress {τ _(xz,s) τ _(yz,s) τ _(zz,s) p _(w) p _(α)} and a displacement {ū_(x) ² ū_(y) ² ū_(z) ² W ² V ²} at any position inside the subgrade in the frequency domain, and a vibration velocity v _(ij)=iωū_(ij) and an acceleration ā_(ij)=ω²ū_(ij) of the subgrade can be calculated according to the displacement;

wherein the processor is further configured to use an inverse fast Fourier transform to convert the solution f of the parameter to be solved in the frequency domain into a solution f in the time domain, wherein the solution f in the time domain is a stress {τ_(xz,s) τ_(yz,s) τ_(zz,s) p_(w) p_(α)}, a displacement {u² _(x) u² _(y) u² _(z) W² V²}, a velocity v_(ij) and an acceleration a_(ij).

In combination of the second aspect of the present disclosure, in the first possible implementation mode according to the second aspect of the present disclosure, the processor is further configured to:

acquire a task of determining the performance of the subgrade of the target expressway section.

In combination of the first possible implementation mode according to the second aspect of the present disclosure, in the second possible implementation mode according to the second aspect of the present disclosure, the processor is specifically configured to:

monitor whether project data on a system are updated; and

in response to monitoring that the project data are updated, acquire the task of determining the performance of the subgrade of the target expressway section based on the updated project data.

In combination of the first possible implementation mode according to the second aspect of the present disclosure, in the third possible implementation mode according to the second aspect of the present disclosure, the processor is specifically configured to:

receive the task of determining the performance of the subgrade of the target expressway section as inputted by a user.

In combination of the second aspect of the present disclosure, in the fourth possible implementation mode according to the second aspect of the present disclosure, the processor further configured to:

output the solution f in the time domain.

In combination of the fourth possible implementation mode according to the second aspect of the present disclosure, in the fifth possible implementation mode according to the second aspect of the present disclosure, the processor is specifically configured to:

output a performance determination report of the subgrade of the target expressway section, wherein the performance determination report includes the solution f in the time domain.

In combination of the second aspect of the present disclosure, in the sixth possible implementation mode according to the second aspect of the present disclosure, the processor is further configured to:

determine whether a value of the solution f in the time domain reaches a predetermined alarm value; and

in response to determining the value of the solution f in the time domain reaches the predetermined alarm value, output an alarm prompt.

In a third aspect, the present disclosure provides a computer readable storage medium which stores a plurality of instructions that are suitable to be loaded by the processor to execute the method provided in the first aspect of the present disclosure or in any possible implementation mode according to the first aspect of the present disclosure.

From the above contents, it can be known that the present disclosure has the following advantageous effects:

It can be seen from the above contents that the present disclosure analyzes the performance of the subgrade of an expressway from the construction of the structural model of the expressway to the derivation of the dynamic governing equations, then removes time terms in the dynamic governing equations and constructs the corresponding dynamic stiffness matrices through the double Fourier transform, and then introduces boundary conditions to continue to solve the dynamic stiffness matrices, and finally converts the parameter solution results in the frequency domain into the parameter solution results in the time domain through the inverse fast Fourier transform, in such a way to accurately calculate the stress, displacement, velocity and acceleration at any position inside the subgrade while analyzing the performance of the subgrade of the target expressway section, in full consideration of the moving speed and vibration characteristics of the vehicle load, the layered characteristics of the expressway and the unsaturated characteristics of the subgrade soil, thereby providing a strong data support for the development of the project.

BRIEF DESCRIPTION OF THE DRAWINGS

In order to illustrate the technical solutions of the embodiments of the present disclosure more clearly, the drawings used in the description of the embodiments will be briefly introduced. Obviously, the drawings in the following description are merely some embodiments of the present disclosure. For those skilled in the art, other drawings can also be obtained from these drawings without creative labor.

FIG. 1 shows a flow chart schematic view of a method for determining a performance of an unsaturated subgrade of a multi-layered expressway section of the present disclosure;

FIG. 2 shows a structural schematic view of an expressway model of the present disclosure;

FIG. 3 shows a structural schematic view of an apparatus for determining a performance of an unsaturated subgrade of a multi-layered expressway section of the present disclosure; and

FIG. 4 shows a structural schematic view of a processing device of the present disclosure.

DETAILED DESCRIPTION

The technical solutions in the embodiments of the present disclosure will be clearly and completely described below with reference to the drawings in the embodiments of the present disclosure. Obviously, the described embodiments are only a part, but not all of the embodiments of the present disclosure. Based on the embodiments in the present disclosure, all other embodiments obtained by those skilled in the art without creative work shall fall within the scope of protection of the present disclosure.

The terms “first”, “second”, etc. in the description and claims, as well as the drawings, of the present disclosure are used to distinguish similar objects, and not necessarily used to describe a specific order or sequence. It should be understood that the data so used can be interchanged where appropriate, so that the embodiments described herein can be implemented in an order other than what is illustrated or described herein. In addition, the terms “comprise” and “have” and any variations thereof are intended to cover non-exclusive inclusion. For example, a process, method, system, product or device that comprises a series of steps or modules are not limited to those steps or modules that are clearly listed, but can include other steps or modules that are not clearly listed or are inherent to these processes, methods, products or devices. The naming or numbering of steps in the present disclosure does not mean that the steps in the process flow must be executed according to the time/logical order indicated by the naming or numbering. The execution order of the named or numbered process steps can be changed according to the technical purpose to be achieved, as long as the same or similar technical effects can be achieved.

The division of modules in the present disclosure is a logical division. There can be another division manner in practical application, for example, a plurality of modules can be combined or integrated in another system, or some features can be ignored or not implemented. In addition, the mutual coupling or direct coupling or communication connection between the displayed or discussed modules can be realized through some interfaces, or the indirect coupling or communication connection between modules can be electrical or other similar forms, which are not limited in the present disclosure. In addition, the modules or sub-modules described as separate components can be or can not be physically separated, can be or can not be physical modules, or can be distributed to a plurality of circuit modules, and some or all of them can be selected according to actual needs to achieve the purpose of the solution of the present disclosure.

Before introducing the method for determining a performance of an unsaturated subgrade of a multi-layered expressway section provided by the present disclosure, the background art involved in the present disclosure shall be introduced first.

The method and the apparatus for determining a performance of an unsaturated subgrade of a multi-layered expressway section, as well as the computer readable storage medium, provided by the present disclosure are applicable to the processing device, which are used to accurately calculate the stress, displacement, velocity and acceleration at any position inside the subgrade while analyzing the performance of the subgrade of the target expressway section, in full consideration of the moving speed and vibration characteristics of the vehicle load, the layered characteristics of the expressway and the unsaturated characteristics of the subgrade soil.

The subject for executing the method for determining a performance of an unsaturated subgrade of a multi-layered expressway section mentioned in the present disclosure can be an apparatus for determining the performance of the subgrade of the expressway section, or different types of processing devices such as servers, physical hosts or user equipment (UE) integrating the apparatus for determining the performance of the subgrade of the expressway section, wherein the apparatus for determining the performance of the subgrade of the expressway section can be implemented in the form of hardware or software, UE can specifically be a terminal device such as a smart phone, tablet, laptop, desktop computer or personal digital assistant (PDA), and the processing device can be arranged in the form of device cluster.

The method for determining a performance of an unsaturated subgrade of a multi-layered expressway section provided by the present disclosure will be introduced.

First, with reference to FIG. 1 , it shows a flow chart schematic view of a method for determining a performance of an unsaturated subgrade of a multi-layered expressway section of the present disclosure. The method for determining the performance of the subgrade of the expressway section provided by the present disclosure can specifically include the following steps.

In step S101: an Odemark equivalent hypothesis of modulus-thickness is used by a processing device to simplify a road structure of a target expressway section into a three-layered model, wherein the three-layered model is

$\begin{matrix} {H_{e} = {\left( \frac{E}{E_{1}} \right)^{1/3}H}} & (1) \end{matrix}$

wherein H_(e) is an equivalent thickness, E is a modulus of a layer to be transformed, H is a thickness of the layer to be transformed, and E₁ is a modulus of a target layer.

It can be understood that reference can be made to a structural schematic view of an expressway model of the present disclosure as shown in FIG. 2 . In the present disclosure, as compared with the prior art which simply simplifies an expressway into a layered elastic medium, the present disclosure will select a flexible pavement or a rigid pavement based on the design scheme or actual situation of the expressway to construct a fine road structural model with the assumption that the base layer is an elastic medium and the subgrade is an unsaturated porous medium, which provides a solid foundation for subsequent detailed performance analysis.

In step S102: based on the three-layered model, a base layer of the target expressway section is set as an elastic medium and the subgrade of the target expressway section is set as an unsaturated porous medium, by the processing device, to derive and obtain a rigid pavement governing equation, a flexible pavement and base layer governing equation, and an unsaturated subgrade governing equation,

wherein the rigid pavement governing equation is:

$\begin{matrix} {{{D_{p}\left\lbrack {\frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial x^{4}} + {2\frac{\partial^{4}{W\left( {x,y,t} \right)}}{{\partial x^{2}}{\partial y^{2}}}} + \frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial y^{4}}} \right\rbrack} + {\rho_{p}h_{p}{\overset{¨}{W}\left( {x,y,t} \right)}} + {\varphi\left( {x,y,t} \right)}} = {f\left( {x,y,t} \right)}} & (2) \end{matrix}$

wherein D_(p)=E_(p)h_(p) ³/[12(1−μ_(p) ²)] is a bending rigidity of a pavement slab, W(x, y, t) is a vertical displacement of the pavement slab, E_(p), ρ_(p) and h_(p) are respectively a modulus, a density and a thickness of the pavement slab, f(x, y, t) is a traffic load on the pavement slab, and φ(x, y, t) is a reaction force exerted by the base layer on the pavement slab.

It can be understood that in the present disclosure, the function (⋅) represents the derivation with respect to time, and also corresponds to the relevant processing of the time terms in the subsequent contents.

The flexible pavement and base layer governing equation is

$\begin{matrix} \left\{ \begin{matrix} {{{G_{b}u_{i,{jj}}^{b}} + {\left( {\lambda_{b} + G_{b}} \right)u_{j,{ji}}^{b}}} = {\rho_{b}{\overset{¨}{u}}_{t}^{b}}} \\ {\tau_{ij}^{b} = {{\lambda_{b}\theta_{b}\delta_{ij}} + {2G_{b}\varepsilon_{ij}^{b}}}} \end{matrix} \right. & (3) \end{matrix}$

wherein u^(b) is a displacement of an elastic medium layer, G_(b) and λ_(b) are Lame constants of the elastic medium layer, ρ_(b) is a density of the elastic medium layer, τ^(b) is a stress of the elastic medium layer, τ^(b) and e^(b) respectively represent the stress of the elastic medium layer and a strain of the elastic medium layer, and δ_(ij) is a Kronecker symbol,

the unsaturated subgrade governing equation is

$\begin{matrix} \left\{ \begin{matrix} {{- p_{w,t}} = {{\rho_{w}{\overset{¨}{u}}_{i}} + {\frac{\rho_{w}}{nS_{r}}{\overset{¨}{W}}_{t}} + {\frac{\rho_{w}g}{k_{w}}{\overset{.}{W}}_{i}}}} \\ {{- p_{a,i}} = {{\rho_{a}{\overset{¨}{u}}_{i}} + {\frac{\rho_{a}}{n\left( {1 - S_{r}} \right)}{\overset{¨}{V}}_{i}} + {\frac{\rho_{a}g}{k_{a}}{\overset{.}{V}}_{i}}}} \\ {{{A_{11}{\overset{.}{p}}_{w}} + {A_{12}{\overset{.}{p}}_{a}} + {A_{13}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{14}{\nabla \cdot {\overset{.}{V}}_{i}}}} = 0} \\ {{{A_{21}{\overset{.}{p}}_{w}} + {A_{22}{\overset{.}{p}}_{a}} + {A_{23}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{24}{\nabla \cdot {\overset{.}{W}}_{i}}}} = 0} \\ {{{\mu{\nabla^{2}{\overset{.}{u}}_{t}}} + {\left( {\lambda + G} \right){\nabla e}} - {a\chi{\nabla p_{w}}} - {{a\left( {1 - \chi} \right)}{\nabla p_{a}}}} = {{\rho{\overset{¨}{u}}_{t}} + {\rho_{w}{\overset{¨}{W}}_{t}} + {\rho_{a}{\overset{¨}{V}}_{t}}}} \end{matrix} \right. & (4) \end{matrix}$

wherein

${A_{11} = {{nA_{s}} + {\left( {1 - S_{r}} \right)\frac{{a\chi} - {nS_{r}}}{K_{s}}}}},$ ${A_{12} = {{\left( {1 - S_{r}} \right)\left\lbrack {\frac{{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}}{K_{s}} + \frac{n}{K_{a}}} \right\rbrack} - {nA_{s}}}},$ ${A_{13} = {\left( {1 - S_{r}} \right)\left( {1 - \frac{K_{b}}{K_{s}}} \right)}},$ A₁₄ = 1, ${A_{21} = {\frac{\left( {{a\chi} - {nS_{r}}} \right)S_{r}}{K_{s}} - {nA_{s}} + \frac{nS_{r}}{K_{w}}}},$ ${A_{22} = {\frac{\left( {\left\lbrack {{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}} \right\rbrack S_{r}} \right.}{\left( K_{s} \right)} + {nA_{s}}}},$ ${A_{23} = {\left( {1 - \frac{K_{b}}{K_{s}}} \right)S_{r}}},$ A₂₄ = 1, ${A_{s} = {{- \alpha}{{md}\left( {1 - {S_{w}0}} \right)}{\left( S_{e} \right)^{\frac{m + 1}{m}}\left\lbrack {\left( S_{e} \right)^{- \frac{1}{m}} - 1} \right\rbrack}^{\frac{d - 1}{d}}}},$ a = 1 − K_(b)/K_(s), K_(b)andK_(s) are bulk compressibility moduli of soil skeletons and soil particles respectively, p_(w) and p_(a) are pressures of pore water and gas respectively, ρ_(w) and ρ_(a) are densities of pore water and gas respectively, W_(i) and V_(i) are displacement components of pore water and pore gas relative to soil particles in the i direction, g represents a gravitational acceleration, k_(w) and k_(a) represent permeability coefficients of pore water and pore gas respectively, x is an effective stress parameter, α, m and d are fitting parameters of a soil-water characteristic curve model, S_(e) is an effective saturation, S_(w0) is an irreducible saturation, n is a soil porosity, and S_(r) is a saturation.

It should be noted herein that in the governing equation (4), the present disclosure introduces an important input index, namely the soil saturation S_(r), and deems that the change of the saturation can have a significant impact on the dynamic response of the overall structure. Thus, the introduction of the saturation S_(r) can facilitate the analysis and processing of the entire performance, and draw attention to the impact of the change of water-containing state inside the subgrade on the overall structure.

Furthermore, as for the governing equation (4), it can be seen that the equation relates to parameters of three aspects, i.e., solid, liquid and gas inside soil, thereby drawing attention to the effect of the coupling between solid, liquid and gas phases in soil, and to the multiphase of soil which is missed in the prior art.

In step S103: based on an action form of a vehicle load, time terms of variables in the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation are respectively set as a simple harmonic form e^(iΩ) ⁰ ^(t) and separated them from the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation by the processing device, and then a double Fourier transform is used by the processing device to respectively reduce partial derivatives in a three-dimensional space domain of the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as ordinary differential equations, to obtain element stiffness matrices of different layer positions after solving and then establish a total dynamic stiffness matrix, wherein the double Fourier transform is

$\begin{matrix} {{{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)} = {\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{f\left( {x,\ y,z,t} \right)e^{- {i({{\beta x} + {\gamma y}})}}{dxdy}}}}},} \\ {{{f\left( {x,y,z,t} \right)} = {\frac{1}{4\pi^{2}}{\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)}e^{i({{\beta x} + {\gamma y}})}d\beta d\gamma}}}}},} \end{matrix}$

wherein β and γ are wave numbers in x and y directions respectively, and f and f are variables in the space domain and corresponding transform domain respectively,

the rigid pavement governing equation in the transform domain is: W (β,γ)[D _(p)(β²+γ²)²−Ω²ρ_(p) h _(p)]+φ(β,γ)= f (β,γ)  (5)

a stress-strain relationship of the flexible pavement in the transform domain is:

$\begin{matrix} {\left\{ {\overset{¯}{\sum}}_{P} \right\} = {{{\left\lbrack S^{p} \right\rbrack\left\lbrack D^{p} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{P} \right\}} = {\left\lbrack K^{p} \right\rbrack\left\{ {\overset{¯}{U}}_{P} \right\}}}} & (6) \end{matrix}$ ${\left\lbrack D^{P} \right\rbrack = \begin{bmatrix} \frac{{- i}\beta e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\beta}{k_{1}^{2}} & e^{{- \alpha_{2}}h} & 1 & 0 & 0 \\ \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\gamma}{k_{1}^{2}} & 0 & 0 & e^{{- \alpha_{2}}h} & 1 \\ \frac{{- i}\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{i\alpha_{1}}{k_{1}^{2}} & {\frac{\beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \beta}{\alpha_{2}} & {\frac{\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \gamma}{\alpha_{2}} \\ \frac{{- i}\beta}{k_{1}^{2}} & \frac{{- i}\beta e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 1 & e^{{- \alpha_{2}}h} & 0 & 0 \\ \frac{{- i}\gamma}{k_{1}^{2}} & \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 0 & 0 & 1 & e^{{- \alpha_{2}}h} \\ \frac{{- i}\alpha_{1}}{k_{1}^{2}} & \frac{i\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{\beta}{\alpha_{2}} & {\frac{- \beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\gamma}{\alpha_{2}} & {\frac{- \gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \end{bmatrix}},$ ${\left\lbrack S^{P} \right\rbrack = {G_{p}\begin{bmatrix} {\frac{2i\beta\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i{\beta\alpha}_{1}}{k_{1}^{2}}} & {{- \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} \\ {\frac{2i{\gamma\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i\gamma\alpha_{1}}{k_{1}^{2}}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} \\ {\frac{i\left( {{2\alpha^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{i\left( {{2\alpha^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}} & {{- 2}\beta e^{{- \alpha_{2}}h}} & {{- 2}\beta} & {{- 2}{\gamma e}^{{- \alpha_{2}}h}} & {{- 2}\gamma} \\ \frac{{- 2}i\beta\alpha_{1}}{k_{1}^{2}} & {\frac{2i{\beta\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} \\ \frac{{- 2}i\gamma\alpha_{1}}{k_{1}^{2}} & {\frac{2i\gamma\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}} & {\frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {2\beta} & {2\beta e^{{- \alpha_{2}}h}} & {2\gamma} & {2\gamma e^{{- \alpha_{2}}h}} \end{bmatrix}}},$

wherein superscript 0 of elements in the matrix represents a top surface of a flexible pavement layer, and superscript 1 of elements in the matrix represents a bottom surface of the flexible pavement layer, a stress vector {Σ _(P)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ τ _(xz,p) ¹ −τ _(yz,p) ¹ −iτ _(zz,p) ¹}^(T), [D^(p)] is a displacement matrix of the flexible pavement, [S^(p)] is a stress matrix of the flexible pavement, a displacement vector {Ū_(P)}={ū_(x,p) ⁰ ū_(y,p) ⁰ iū_(z,p) ⁰ ū_(x,p) ¹ ū_(y,p) ¹ iū_(z,p) ¹}^(T), [K^(p)]=[S^(p)][D^(p)]⁻¹ s a dynamic stiffness matrix of the flexible pavement layer, α₁ ²=β²+γ²−k₁ ², α₂ ²=β²+γ²−k₂ ², k₁ ²=ω²/c₁ ², k₂ ²=ω²/c₂ ², c is a traffic load velocity, c₁ and c₂ are a compression wave velocity and a shear wave velocity of the elastic medium respectively, ω=Ω₀−βc, h represents a thickness of the flexible pavement layer, and i is an imaginary number,

a stress-strain relationship of the elastic base layer in the transform domain is: {Σ _(b) }=[S ^(b) ][D ^(b)]⁻¹ {Ū _(b) }=[K ^(b) ]{Ū _(b)}  (7)

a stress-strain relationship of the unsaturated subgrade in the transform domain is:

$\begin{matrix} {\left\{ {\overset{¯}{\sum}}_{s} \right\} = {{{\left\lbrack S^{s} \right\rbrack\left\lbrack D^{s} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{s} \right\}} = {\left\lbrack K^{s} \right\rbrack\left\{ {\overset{¯}{U}}_{s} \right\}}}} & (8) \end{matrix}$ ${\left\lbrack D^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma}{\beta}}e^{{- \lambda_{0}}z}} & {\frac{i\lambda_{0}}{\beta}e^{{- \lambda_{0}}z}} & {S_{1}i\beta e^{{- \lambda_{1}}z}} & {S_{2}i\beta e^{{- \lambda_{2}}z}} & {S_{3}i\beta e^{{- \lambda_{3}}z}} \\ e^{{- \lambda_{0}}z} & 0 & {S_{1}i\gamma e^{{- \lambda_{1}}z}} & {S_{2}i\gamma e^{{- \lambda_{2}}z}} & {S_{3}i\gamma e^{{- \lambda_{2}}z}} \\ 0 & {- {ie}^{{- \lambda_{0}}z}} & {{- i}S_{1}\lambda_{1}e^{{- \lambda_{1}}z}} & {{- i}S_{2}\lambda_{2}e^{{- \lambda_{1}}z}} & {{- i}S_{3}\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{w}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{w1}}{b^{w}\rho_{w}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{w2}}{b^{w}\rho_{w}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{w3}}{b^{w}\rho_{w}}}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{a}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{a1}}{b^{a}\rho_{a}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{a2}}{b^{a}\rho_{a}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{a3}}{b^{a}\rho_{a}}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$

${\left\lbrack S^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma\lambda_{0}G}{\beta}}e^{{- \lambda_{0}}z}} & {{{iG}\left( {\frac{\lambda_{0}^{2}}{\beta} + \beta} \right)}e^{{- \lambda_{0}}z}} & {2GS_{1}i\beta\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i\beta\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i\beta\lambda_{3}e^{{- \lambda_{3}}z}} \\ {G\lambda_{0}e^{{- \lambda_{0}}z}} & {iG\gamma e^{{- \lambda_{0}}z}} & {2GS_{1}i{\gamma\lambda}_{1}e^{{- \lambda_{2}}z}} & {2GS_{2}i\gamma\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i\gamma\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- i}2G\lambda_{0}e^{{- \lambda_{0}}z}} & {{- i}B_{1}e^{{- \lambda_{1}}z}} & {{- i}B_{2}e^{{- \lambda_{2}}z}} & {{- i}B_{3}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{w1}}e^{{- \lambda_{1}}z}} & {{- f_{w2}}e^{\lambda_{2}z}} & {{- f_{w3}}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{a1}}e^{{- \lambda_{1}}z}} & {{- f_{a2}}e^{{- \lambda_{2}}z}} & {{- f_{a3}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$

wherein

${\left\{ {\sum\limits^{\_}}_{s} \right\} = \begin{Bmatrix} {- {\overset{¯}{\tau}}_{{xz},s}\ } & {- {\overset{¯}{\tau}}_{{yz},s}} & {{- i}{\overset{¯}{\tau}}_{{zz},s}} & {- {\overset{¯}{p}}_{w}\ } & {- {\overset{¯}{p}}_{a}} \end{Bmatrix}^{T}},$ ${\left\{ {\overset{¯}{U}}_{s} \right\} = \begin{Bmatrix} {\overset{¯}{u}}_{x,s} & {\overset{¯}{u}}_{y,s} & {i{\overset{¯}{u}}_{z,s}} & \overset{¯}{W} & \overset{¯}{V} \end{Bmatrix}^{T}},$ B_(wn) = ω²ρ_(w)S_(n)λ_(n) − f_(wn)λ_(n), B_(an) = ω²ρ_(a)S_(n)λ_(n) − f_(an)λ_(n), B_(n) = 2μS_(n)λ_(n)² + λ − aχf_(wn) − a(1 − χ)f_(an), ${S_{n} = {- \frac{\lambda + \mu + {b_{1}f_{wn}} + {b_{2}f_{an}}}{b_{3} + {\mu d_{n}}}}},$ ${f_{an} = \frac{\left( {{b_{13}b_{31}} - {b_{11}b_{33}} + {b_{11}d_{n}}} \right)}{\left( {{b_{11}b_{32}} - {b_{31}b_{12}} + {b_{31}d_{n}}} \right)}},$ ${f_{wn} = \frac{\left( {{b_{32}b_{23}} - {b_{22}b_{33}} + {b_{22}d_{n}}} \right)}{\left( {{b_{22}b_{31}} - {b_{32}b_{21}} + {b_{32}d_{n}}} \right)}},$ ${\lambda_{n} = \sqrt{\beta^{2} + \gamma^{2} + d_{n}}},d_{n}$ is roots of the equation

d_(n)³ + b₄d_(n)² + b₅d_(n) + b₆ = 0, b₄ = −b₁₂ − b₃₃ − b₂₁, b₅ = b₁₂b₃₃ − b₁₃b₃₂ + b₂₁b₁₂ + b₂₁b₃₃ − b₂₂b₁₁ − b₂₃b₃₁ b₆ = −b₂₁b₁₂b₃₃ − b₁₃b₃₁b₂₂ − b₂₃b₁₁b₃₂ + b₂₁b₁₃b₃₂ + b₂₂b₁₁b₃₃ + b₂₃b₃₁b₁₂, $\begin{matrix} {{n = 1},2,3,} & {{b_{1} = {{- {\alpha\chi}} - \frac{\omega^{2}}{b^{w}}}},} \end{matrix}$ ${b_{2} = {\left( {1 - \chi} \right) - \frac{\omega^{2}}{b^{a}}}},$ ${b_{3} = {{\rho\omega}^{2} + \frac{\omega^{4}\rho_{w}}{b^{w}} + \frac{\omega^{4}\rho_{a}}{b^{a}}}},$ b₁₁ = b^(a)ρ_(a)A₁₁, b₁₂ = b^(a)ρ_(a)A₁₂, b₁₃ = ω²ρ_(a) + b^(a)ρ_(a)A₁₃, b₂₁ = b^(w)ρ_(w)A₂₁, b₂₂ = b^(w)ρ_(w)A₂₂, b₂₃ = ω²ρ_(w) + b^(w)ρ_(w)A₂₃, ${b_{31} = {- \frac{{b_{1}b_{21}} + {b_{2}b_{11}}}{\lambda + {2\mu}}}},$ ${b_{32} = {- \frac{{b_{1}b_{22}} + {b_{2}b_{12}}}{\lambda + {2\mu}}}},$ ${b_{33} = {- \frac{{b_{1}b_{23}} + {b_{2}b_{13}} + b_{3}}{\lambda + {2\mu}}}},$ ${b^{w} = {{- \frac{\omega^{2}}{nS_{r}}} + \frac{i{\omega g}}{k_{w}}}},$ ${b^{a} = {{- \frac{\omega^{2}}{n\left( {1 - S_{r}} \right)}} + \frac{i\omega g}{k_{a}}}},$

a total stiffness matrix of a flexible pavement expressway model is: {Σ _(f) }=[K ^(f) ]{Ū _(f)}  (9)

wherein {Σ _(f)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ 0 0 0 0 0 0 0 0}^(T), {Ū_(f)}={ū_(x) ⁰ ū_(y) ⁰ iū_(z) ⁰ ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), superscripts 0, 1, 2 correspond to different layer positions shown from top to bottom in the three-layered model (the model structure as shown in FIG. 2 ), and [K^(f)] is a matrix of 11×11,

a total stiffness matrix of a rigid pavement expressway model is: {Σ _(r) }=[K ^(r) ]{Ū _(r)}  (10)

wherein {Σ _(r)}={−τ _(xz,b) ¹ −τ _(yz,b) ¹ −iτ _(zz,b) ¹ 0 0 0 0 0}^(T), {Ū_(r)}={ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), and [K^(r)] is a matrix of 8×8;

It can be understood that the construction of the stiffness matrix is used to provide the analysis unit of finite elements for subsequent performance analysis.

In step S104: a load form and a load size are determined by the processing device according to a specific vehicle, and the load form and the load size are substituted as boundary conditions into the stiffness matrix of a corresponding model by the processing device to solve and obtain a solution f of a parameter to be solved in a frequency domain, wherein the load form and the load size determined according to the specific vehicle are

${f\left( {x_{1},y_{1},t} \right)} = \left\{ {\begin{matrix} {\frac{{Pe}^{i\Omega_{0}t}}{40l_{1}l_{2}}:\left\{ \begin{matrix} {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + l_{2}} \right)}} \end{matrix} \right.} \\ {0;{others}} \end{matrix},} \right.$

wherein P represents a total axle load of the specific vehicle, and in the moving coordinate system x=x₁−Vt, y=y₁, z=z₁.

In the transform domain, the load form and the load size determined according to the specific vehicle are

$\begin{matrix} {{{\overset{¯}{f}\left( {\beta,\gamma} \right)} = \frac{P\left( {e^{{- \gamma}d_{w}} + e^{i\gamma d_{w}}} \right)}{10l_{1}l_{2}\beta\gamma}}\text{ }\left\lbrack {{{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}e^{{- i}\beta d_{c}}} + {{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}\left( {e^{{- i}\beta d_{a}} + e^{{- i}\beta{d}_{a}}} \right)\left( {e^{{- i}\gamma d_{t}} + e^{i\gamma d_{t}}} \right)}} \right\rbrack} & (12) \end{matrix}$

the boundary condition of the rigid pavement expressway model is:

$\begin{matrix} \left\{ \begin{matrix} {{\tau_{{zz},b}^{1}\left( {x,y,t} \right)} = {- {\varphi\left( {x,y,t} \right)}}} \\ {{\tau_{{xz},b}^{1}\left( {x,y,t} \right)} = {{\tau_{{yz},b}^{1}\left( {x,\ y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \\ {{u_{{zz},b}^{1}\left( {x,y,t} \right)} = {W\left( {x,y,t} \right)}} \end{matrix} \right. & (13) \end{matrix}$

the boundary condition of the rigid pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\begin{matrix} \left\{ \begin{matrix} {{{\overset{\_}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {- {\overset{\_}{\varphi}\left( {\beta,\gamma} \right)}}} \\ {{{\overset{\_}{\tau}}_{{xz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{\tau}}_{{yz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \\ {{{\overset{\_}{u}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {W\left( {\beta,\gamma} \right)}} \end{matrix} \right. & (14) \end{matrix}$

wherein

${{\overset{¯}{\varphi}\left( {\beta,\gamma} \right)} = \frac{\overset{¯}{f}\left( {\beta,\gamma} \right)}{1 + {B_{33}{D_{P}\left( {\beta^{2} + \gamma^{2}} \right)}^{2}} - {B_{33}\omega^{2}m_{b}}}},$ and B₃₃ is an element in a third row and a third column of [K^(r)]⁻¹ and corresponds to a vertical displacement of the base layer,

the boundary condition of the flexible pavement expressway model is:

$\begin{matrix} \left\{ \begin{matrix} {{\tau_{{zz},p}^{0}\left( {x,y,t} \right)} = {f\left( {x,y,t} \right)}} \\ {{\tau_{xz}^{0}\left( {x,y,t} \right)} = {{\tau_{yz}^{0}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \end{matrix} \right. & (15) \end{matrix}$

the boundary condition of the flexible pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\begin{matrix} \left\{ \begin{matrix} {{{\overset{\_}{\tau}}_{{zz},p}^{0}\left( {\beta,\gamma} \right)} = {f\left( {\beta,\gamma} \right)}} \\ {{{\overset{\_}{\tau}}_{xz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{\tau}}_{yz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \end{matrix} \right. & (16) \end{matrix}$

wherein for the rigid pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the rigid pavement expressway model, and for the flexible pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the flexible pavement expressway model, in such a way to obtain the displacement of each surface, the displacement {ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²} of the top surface of the subgrade is substituted into the stress-strain relationship of the unsaturated subgrade in the transform domain to obtain a stress {τ _(xz,s) τ _(yz,s) τ _(zz,s) p _(w) p _(α)} and a displacement {ū_(x) ² ū_(y) ² ū_(z) ² W ² V ²} at any position inside the subgrade in the frequency domain, and a vibration velocity and an acceleration of the subgrade can be calculated according to the displacement, the vibration velocity and the acceleration are v _(ij)=iωū_(ij) and ā_(ij)=−ωū_(ij).

In this process, the load of the vehicle that may be involved in the actual expressway environment, such as a heavy loaded 10-wheeled truck, is introduced into the stiffness matrix constructed in the Step S104 for solving. During the solving process, the constraint of boundary conditions of the configuration may be involved. The solution range is constrained so that the solution results may be more accurate.

It should be understood that the solution results at this stage are still parameters in the frequency domain, and further data processing is still required to facilitate practical use.

In Step S105: an inverse fast Fourier transform is used by the processing device to convert the solution f of the parameter to be solved in the frequency domain into a solution f in the time domain, wherein the solution f in the time domain is a stress {τ_(xz,s) τ_(yz,s) τ_(zz,s) p_(w) p_(α)}, a displacement {u² _(x) u² _(y) u² _(z) W² V²}, a velocity v_(ij) and an acceleration a_(ij).

Here, corresponding to the above-mentioned double Fourier transform, parameters in the frequency domain solved in Step S104 are converted into parameters in the time domain by the inverse fast Fourier transform, which are the final required parameters and also the parameters that can be used conveniently in practical applications.

For instance, 2048×2048 points can be selected between −16<β<16 m−1 and −16<γ<16 m−1, and values of parameters to be solved in the frequency domain at each integration point can be calculated.

It can be seen from the above contents that the present disclosure analyzes the performance of the subgrade of the expressway from the construction of the structural model of the expressway to the derivation of the dynamic governing equations, then removes time terms in the dynamic governing equations and constructs corresponding dynamic stiffness matrices through double Fourier transform, and then introduces boundary conditions to continue to solve the dynamic stiffness matrices, and finally converts the parameter solution results in the frequency domain into the parameter solution results in the time domain through the inverse fast Fourier transform, in such a way to accurately calculate the stress, displacement, velocity and acceleration at any position inside the subgrade while analyzing the performance of the subgrade of the target expressway section, in full consideration of the moving speed and vibration characteristics of the vehicle load, the layered characteristics of the expressway and the unsaturated characteristics of the subgrade soil, thereby providing a strong data support for the development of the project.

Furthermore, each step of the above embodiment as shown in FIG. 1 and possible implementation modes thereof in practical applications will be explained in detail in combination with practical applications.

It can be understood that as the trigger of the analysis and processing of the performance of the subgrade of the expressway section, it can also be realized in the form of a task in practical applications.

Namely, the method for determining the performance of the subgrade of the expressway section provided by the present disclosure further includes Step S106, which is prior to Step S101:

the processing device acquiring a task of determining the performance of the subgrade of the target expressway section.

Wherein the acquisition of the task of determining the performance of the subgrade of the target expressway section can be actively acquired or passively received, which may be adjusted as actually required.

As an actively acquired implementation mode, Step S106 specifically includes:

the processing device monitoring whether project data on a system are updated; and

in response to monitoring that the project data are updated, the processing device acquiring the task of determining the performance of the subgrade of the target expressway section according to the updated target project data.

It can be understood that the processing device itself can be a device included in the office system of the project for the purpose of monitoring and updating within the system.

Alternatively, the processing device may also externally establish a communication connection with the office system of the project so as to interact with the system for the purpose of monitoring and updating.

As a passively received implementation mode, Step S106 specifically comprises:

the processing device receiving the task of determining the performance of the subgrade of the target expressway section as inputted by a user.

Wherein the task can be specifically the one inputted into the device locally by a user of human-computer interaction, which is applicable to on-site office scenario.

Alternatively, the task can also be the one inputted by a user of human-computer interaction and received by the device locally from another device, which is applicable to remote task assigning scenario.

Furthermore, after the parameter solution results in the time domain are obtained through conversion in Step S105, they can also be outputted, namely, the method for determining the performance of the subgrade of the expressway section can further include Step S107:

the processing device outputting the solution f in the time domain.

It can be understood that the output processing herein is executed under a predetermined output strategy. The output carrier involved can be hardware such as a display screen, or even a loudspeaker, or an indicator light, and the output manner can be specifically a text message, instant message, document, image, voice, vibration, lighting effect, etc., which can be configured according to the output effect demands in practical applications.

The present disclosure also provides a practical implementation mode in combination with the project scenario in practical applications. In the output process, Step S107 can specifically include:

the processing device outputting a performance determination report of the subgrade of the target expressway section, wherein the performance determination report includes the solution f in the time domain.

It can be understood that the solution f in the time domain (the parameter solution result in the time domain) is directly carried in the performance determination report in the original application scenario, which is greatly convenient for application. The report also contains the contents in other different aspects.

Of course, the performance determination report can also be the one arranged for the report working mode that commonly exists in the original application scenario, and can be conveniently embedded into the specific application scenario in a terminal data output processing step.

In addition, for the parameter solution result in the time domain, an early warning mechanism may also be introduced in practical applications. Therefore, the method for determining the performance of the subgrade of the expressway section can also include Step S108, which is subsequent to Step S105:

the processing device determining whether a value of the solution f in the time domain reaches a predetermined alarm value; and

in response to determining the value of the solution f in the time domain reaches the predetermined alarm value, the processing device outputting an alarm prompt.

It can be understood that the predetermined alarm value can be a numerical threshold value, or a time element can also be introduced as the alarm value for dynamic monitoring, such as the change rate in unit time or the change range in unit time, so that different specific warning conditions can be set according to different warning demands.

Then, alarm prompts can also be outputted in different output forms according to actual needs. The output carrier involved can be hardware such as a display screen, or even a loudspeaker, or an indicator light, and the output manner can be specifically a text message, instant message, document, image, voice, vibration, lighting effect, etc., which can be configured according to the output effect demands in practical applications.

The above is the introduction of the method for determining the performance of the subgrade of the expressway section provided by the present disclosure. In order to better implement the method for determining the performance of the subgrade of the expressway section provided by the present disclosure, the present disclosure also provides an apparatus for determining a performance of an unsaturated subgrade of a multi-layered expressway section from the perspective of functional modules.

With reference to FIG. 3 , FIG. 3 shows a structural schematic view of an apparatus for determining a performance of an unsaturated subgrade of a multi-layered expressway section of the present disclosure. In the present disclosure, the apparatus 300 for determining a performance of an unsaturated subgrade of a multi-layered expressway section can specifically include a structural model simplifying unit 301, a governing equation deriving unit 302, a stiffness matrix constructing unit 303, a solving unit 304 and a converting unit 305.

The structural model simplifying unit 301 is configured to use an Odemark equivalent hypothesis of modulus-thickness to simplify a road structure of a target expressway section into a three-layered model, wherein the three-layered model is

${H_{e} = {\left( \frac{E}{E_{1}} \right)^{1/3}\mspace{14mu} H}},$ wherein H_(e) is an equivalent thickness, E is a modulus of a layer to be transformed, H is a thickness of the layer to be transformed, and E₁ is a modulus of a target layer.

The governing equation deriving unit 302 is configured to, based on the three-layered model, set a base layer of the target expressway section as an elastic medium and set the subgrade of the target expressway section as an unsaturated porous medium to derive and obtain a rigid pavement governing equation, a flexible pavement and base layer governing equation and an unsaturated subgrade governing equation, wherein the rigid pavement governing equation is:

${{{D_{p}\left\lbrack {\frac{\partial^{2}{W\left( {x,y,t} \right)}}{\partial x^{4}} + {2\frac{\partial^{4}{W\left( {x,y,t} \right)}}{{\partial x^{2}}{\partial y^{2}}}} + \frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial y^{4}}} \right\rbrack} + {\rho_{p}h_{p}{\overset{¨}{W}\left( {x,y,t} \right)}} + {\varphi\left( {x,y,t} \right)}} = {f\left( {x,y,t} \right)}},$

wherein D_(p)=E_(p)h_(p) ³/[12(1−μ_(p) ²)] is a bending rigidity of a pavement slab, W(x, y, t) is a vertical displacement of the pavement slab, E_(p), ρ_(p) and h_(p) are respectively a modulus, a density and a thickness of the pavement slab, f(x, y, t) is a traffic load on the pavement slab, and φ(x, y, t) is a reaction force exerted by the base layer on the pavement slab,

the flexible pavement and base layer governing equation is

$\left\{ {\begin{matrix} {{{G_{b}u_{i,{jj}}^{b}} + {\left( {\lambda_{b} + G_{b}} \right)u_{j,{ji}}^{b}}} = {\rho_{b}{\overset{¨}{u}}_{i}^{b}}} \\ {\tau_{ij}^{b} = {{\lambda_{b}\theta_{b}\delta_{ij}} = {2\; G_{b}ɛ_{ij}^{b}}}} \end{matrix},} \right.$

wherein u^(b) is a displacement of an elastic medium layer, G_(b) and λ_(b) are Lame constants of the elastic medium layer, ρ_(b) is a density of the elastic medium layer, τ^(b) is a stress of the elastic medium layer, τ^(b) and ε^(b) respectively represent the stress of the elastic medium layer and a strain of the elastic medium layer, and δ_(ij) is a Kronecker symbol,

the unsaturated subgrade governing equation is

$\left\{ {\begin{matrix} {{- p_{w,i}} = {{\rho_{w}{\overset{¨}{u}}_{i}} + {\frac{\rho_{w}}{{nS}_{r}}{\overset{¨}{W}}_{l}} + {\frac{\rho_{w}g}{k_{w}}{\overset{.}{W}}_{i}}}} \\ {{- p_{a,i}} = {{\rho_{a}{\overset{¨}{u}}_{i}} + {\frac{\rho_{a}}{n\left( {1 - S} \right)}{\overset{¨}{V}}_{i}} + {\frac{\rho_{a}g}{k_{a}}{\overset{.}{V}}_{i}}}} \\ {{{A_{11}{\overset{.}{p}}_{w}} + {A_{12}{\overset{.}{p}}_{a}} + {A_{13}{\nabla{\cdot {\overset{.}{u}}_{i}}}} + {A_{14}{\nabla{\cdot {\overset{.}{V}}_{i}}}}} = 0} \\ {{{A_{21}{\overset{.}{p}}_{w}} + {A_{22}{\overset{.}{p}}_{a}} + {A_{23}{\nabla{\cdot {\overset{.}{u}}_{i}}}} + {A_{24}{\nabla{\cdot {\overset{.}{W}}_{i}}}}} = 0} \\ {{{\mu\;{\nabla^{2}{\overset{.}{u}}_{i}}} + {\left( {\lambda + G} \right){\nabla\; e}} - {a\;\mathcal{X}\;{\nabla\; p_{w}}} - {{a\left( {1 - \mathcal{X}} \right)}{\nabla\; p_{a}}}} = {{\rho{\overset{¨}{u}}_{i}} + {\rho_{w}{\overset{¨}{W}}_{i}} + {\rho_{a}{\overset{¨}{V}}_{i}}}} \end{matrix},} \right.$

wherein

${A_{11} = {{n\; A_{s}} + {\left( {1 - S_{r}} \right)\frac{{a\;\mathcal{X}} - {n\; S_{r}}}{K_{s}}}}},$

${A_{12} = {{\left( {1 - S_{r}} \right)\left\lbrack {\frac{{a\left( {1 - \mathcal{X}} \right)} - {n\left( {1 - S_{r}} \right)}}{K_{s}} + \frac{n}{K_{a}}} \right\rbrack} - {n\; A_{s}}}},{A_{13} = {\left( {1 - S_{r}} \right)\left( {1 - \frac{K_{b}}{K_{s}}} \right)}},{A_{14} = 1},{A_{21} = {\frac{\left( {{a\;\mathcal{X}} - {n\; S_{r}}} \right)S_{r}}{K_{s}} - {n\; A_{s}} + \frac{n\; S_{r}}{K_{w}}}},{A_{22} = {\frac{\left\lbrack {{a\left( {1 - \mathcal{X}} \right)} - {n\left( {1 - S_{r}} \right)}} \right\rbrack S_{r}}{K_{s}} + {n\; A_{s}}}},{A_{23} = {\left( {1 - \frac{K_{b}}{K_{s}}} \right)S_{r}}},{A_{24} = 1},{A_{s} = {{- \alpha}\mspace{11mu} m\;{d\left( {1 - S_{w\; 0}} \right)}{\left( S_{e} \right)^{\frac{m + 1}{m}}\left\lbrack {\left( S_{e} \right)^{- \frac{1}{m}} - 1} \right\rbrack}^{\frac{d - 1}{d}}}},{a = {1 - {K_{b}/K_{s}}}},{K_{b}\mspace{14mu}{and}\mspace{14mu} K_{s}\mspace{14mu}{are}\mspace{14mu}{bulk}}$ compressibility moduli of soil skeletons and soil particles respectively, p_(w) and p_(a) are pressures of pore water and gas respectively, ρ_(w) and ρ_(a) are densities of pore water and gas respectively, W_(i) and V_(i) are displacement components of pore water and pore gas relative to soil particles in an i direction, g represents a gravitational acceleration, k_(w) and k_(a) represent permeability coefficients of pore water and pore gas respectively, x is an effective stress parameter, α, m and d are fitting parameters of a soil-water characteristic curve model, S_(e) is an effective saturation, S_(w0) is an irreducible saturation, n is a soil porosity, and S_(r) is a saturation.

The stiffness matrix constructing unit 303 is configured to, based on an action form of a vehicle load, respectively set time terms of variables in the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as a simple harmonic form e^(iΩ) ⁰ ^(t) and separate the simple harmonic form e^(iΩ) ⁰ ^(t) of the time terms from the rigid pavement governing equation, the flexible pavement and base layer governing equation, and the unsaturated subgrade governing equation, and then use a double Fourier transform to respectively reduce partial derivatives in a three-dimensional space domain of the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as ordinary differential equations, to obtain element stiffness matrices of different layer positions after solving and then establish a total dynamic stiffness matrix, wherein the double Fourier transform is

${{\overset{\_}{f}\left( {\beta,\gamma,z,t} \right)} = {\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{f\left( {x,y,z,t} \right)}e^{- {i{({{\beta\mspace{11mu} x} + {\gamma\mspace{11mu} y}})}}}{dxdy}}}}},{{f\left( {x,y,z,t} \right)} = {\frac{1}{4\;\pi^{2}}{\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{\overset{\_}{f}\left( {\beta,\gamma,z,t} \right)}e^{i{({{\beta\; x} + {\gamma\; y}})}}d\;\beta\; d\;\gamma}}}}},$

wherein β and γ are wave numbers in x and y directions respectively, and f and f are variables in the space domain and a corresponding transform domain respectively,

the rigid pavement governing equation in the transform domain is: W (β,γ)[D _(p)(β²+γ²)²−Ω²ρ_(p) h _(p)]+φ(β,γ)= f (β,γ),

a stress-strain relationship of the flexible pavement in the transform domain is:

${\left\{ {\sum\limits^{\_}}_{P}\; \right\} = {{{\left\lbrack S^{p} \right\rbrack\left\lbrack D^{p} \right\rbrack}^{- 1}\left\{ {\overset{\_}{U}}_{P} \right\}} = {\left\lbrack K^{p} \right\rbrack\left\{ {\overset{\_}{U}}_{P} \right\}}}},{\left\lbrack D^{P} \right\rbrack = \begin{bmatrix} \frac{{- i}\;\beta\; e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\;\beta}{k_{1}^{2}} & e^{{- \alpha_{2}}h} & 1 & 0 & 0 \\ \frac{{- i}\;\gamma\; e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\;\gamma}{k_{1}^{2}} & 0 & 0 & e^{{- \alpha_{2}}h} & 1 \\ \frac{{- i}\;\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{i\;\alpha_{1}}{k_{1}^{2}} & {\frac{\beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \beta}{\alpha_{2}} & {\frac{\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \gamma}{\alpha_{2}} \\ \frac{{- i}\;\beta}{k_{1}^{2}} & \frac{{- i}\;\beta\; e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 1 & e^{{- \alpha_{2}}h} & 0 & 0 \\ \frac{{- i}\;\gamma}{k_{1}^{2}} & \frac{{- i}\;\gamma\; e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 0 & 0 & 1 & e^{{- \alpha_{2}}h} \\ \frac{{- i}\;\alpha_{1}}{k_{1}^{2}} & \frac{i\;\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{\beta}{\alpha_{2}} & {\frac{- \beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\gamma}{\alpha_{2}} & {\frac{- \gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \end{bmatrix}},$

${\left\lbrack S^{P} \right\rbrack = {G_{p}\begin{bmatrix} {\frac{2\; i\;\beta\;\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2\; i\;\beta\;\alpha_{1}}{k_{1}^{2}}} & {{- \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} \\ {\frac{2\; i\;\gamma\;\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2\; i\;\gamma\;\alpha_{1}}{k_{1}^{2}}} & {{- \frac{\beta\;\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} \\ {\frac{i\left( {{2\;\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{i\left( {{2\;\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}} & {{- 2}\;\beta\; e^{{- \alpha_{2}}h}} & {{- 2}\beta} & {{- 2}\gamma\; e^{{- \alpha_{2}}h}} & {{- 2}\;\gamma} \\ \frac{{- 2}\; i\;\beta\;\alpha_{1}}{k_{1}^{2}} & {\frac{2\; i\;\beta\;\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} \\ \frac{{- 2}\; i\;\gamma\;\alpha_{1}}{k_{1}^{2}} & {\frac{2\; i\;\gamma\;\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\beta\;\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{{- i}\;\left( {{2\;\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}} & {\frac{- {i\left( {{2\;\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {2\;\beta} & {2\;\beta\; e^{{- \alpha_{2}}h}} & {2\;\gamma} & {2\gamma\; e^{{- \alpha_{2}}h}} \end{bmatrix}}},$

wherein superscript 0 of elements in the matrix represents a top surface of a flexible pavement layer, and superscript 1 of elements in the matrix represents a bottom surface of the flexible pavement layer, a stress vector {Σ _(P)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ τ _(xz,p) ¹ −τ _(yz,p) ¹ −iτ _(zz,p) ¹}^(T), [D^(p)] is a displacement matrix of the flexible pavement, [S^(p)] is a stress matrix of the flexible pavement, a displacement vector {Ū_(P)}={ū_(x,p) ⁰ ū_(y,p) ⁰ iū_(z,p) ⁰ ū_(x,p) ¹ ū_(y,p) ¹ iū_(z,p) ¹}^(T), [K^(p)]=[S^(p)][D^(p)]⁻¹ is a dynamic stiffness matrix of the flexible pavement layer, α₁ ²=β²+γ²−k₁ ², α₂ ²=β²+γ²−k₂ ², k₁ ²=ω²/c₁ ², k₂ ²=ω²/c₂ ², c is a traffic load velocity, c₁ and c₂ are a compression wave velocity and a shear wave velocity of the elastic medium respectively, ω=Ω₀−βc, h represents a thickness of the flexible pavement layer, and i is an imaginary number,

a stress-strain relationship of the elastic base layer in the transform domain is: {Σ _(b) }=[S ^(b) ][D ^(b)]⁻¹ {Ū _(b) }=[K ^(b) ]{Ū _(b)},

a stress-strain relationship of the unsaturated subgrade in the transform domain is:

${\left\{ {\sum\limits^{\_}}_{s} \right\} = {{{\left\lbrack S^{s} \right\rbrack\left\lbrack D^{s} \right\rbrack}^{- 1}\left\{ {\overset{\_}{U}}_{s} \right\}} = {\left\lbrack K^{s} \right\rbrack\left\{ {\overset{\_}{U}}_{s} \right\}}}},{\left\lbrack D_{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma}{\beta}}e^{{- \lambda_{0}}z}} & {\frac{i\;\lambda_{0}}{\beta}e^{{- \lambda_{0}}z}} & {S_{1}i\;\beta\; e^{{- \lambda_{1}}z}} & {S_{2}i\;\beta\; e^{{- \lambda_{2}}z}} & {S_{3}i\;\beta\; e^{{- \lambda_{3}}z}} \\ e^{{- \lambda_{0}}z} & 0 & {S_{1}i\;\gamma\; e^{{- \lambda_{1}}z}} & {S_{2}i\;\gamma\; e^{{- \lambda_{2}}z}} & {S_{3}i\;\gamma\; e^{{- \lambda_{3}}z}} \\ 0 & {- {ie}^{{- \lambda_{0}}z}} & {{- {iS}_{1}}\lambda_{1}e^{{- \lambda_{1}}z}} & {{- i}\; S_{2}\lambda_{2}e^{{- \lambda_{2}}z}} & {{- i}\; S_{3}\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{w}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{w\; 1}}{b^{w}\rho_{w}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{w\; 2}}{b^{w}\rho_{w}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{w\; 3}}{b^{w}\rho_{w}}}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{a}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{a\; 1}}{b^{a}\rho_{a}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{a\; 2}}{b^{a}\rho_{a}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{a\; 3}}{b^{a}\rho_{a}}}e^{{- \lambda_{3}}z}} \end{bmatrix}},{\left\lbrack S^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma\;\lambda_{0}G}{\beta}}e^{{- \lambda_{0}}z}} & {{{iG}\left( {\frac{\lambda_{0}^{2}}{\beta} + \beta} \right)}e^{{- \lambda_{0}}z}} & {2\; G\; S_{1}i\;\beta\;\lambda_{1}e^{{- \lambda_{1}}z}} & {2\; G\; S_{2}i\;\beta\;\lambda_{2}e^{{- \lambda_{2}}z}} & {2\; G\; S_{3}i\;{\beta\lambda}_{3}e^{{- \lambda_{3}}z}} \\ {G\;\lambda_{0}e^{{- \lambda_{0}}z}} & {{iG}\;\gamma\; e^{{- \lambda_{0}}z}} & {2\; G\; S_{1}i\;{\gamma\lambda}_{1}e^{{- \lambda_{1}}z}} & {2\; G\; S_{2}i\;{\gamma\lambda}_{2}e^{{- \lambda_{2}}z}} & {2\; G\; S_{3}i\;{\gamma\lambda}_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- i}\; 2\; G\;\lambda_{0}e^{{- \lambda_{0}}z}} & {{- i}\; B_{1}e^{{- \lambda_{1}}z}} & {{- i}\; B_{2}e^{{- \lambda_{2}}z}} & {{- i}\; B_{3}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{w\; 1}}e^{{- \lambda_{1}}z}} & {{- f_{w\; 2}}e^{{- \lambda_{2}}z}} & {{- f_{w\; 3}}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{a\; 1}}e^{{- \lambda_{1}}z}} & {{- f_{a\; 2}}e^{{- \lambda_{2}}z}} & {{- f_{a\; 3}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$

wherein

${\left\{ {\sum\limits^{\_}}_{s} \right\} = \left\{ {{- {\overset{\_}{\tau}}_{{xz},s}}\mspace{31mu} - {\overset{\_}{\tau}}_{{yz},s}\mspace{31mu} - {i\;{\overset{\_}{\tau}}_{{zz},s}}\mspace{31mu} - {\overset{\_}{p}}_{w}\mspace{31mu} - {\overset{\_}{p}}_{a}} \right\}^{T}},{\left\{ {\overset{\_}{U}}_{s} \right\} = \left\{ {{\overset{\_}{u}}_{x,s}\mspace{31mu}{\overset{\_}{u}}_{y,s}\mspace{31mu} i\;{\overset{\_}{u}}_{z,s}\mspace{31mu}\overset{\_}{W}\mspace{31mu}\overset{\_}{V}} \right\}^{T}},{B_{wn} = {{\omega^{2}\rho_{w}S_{n}\lambda_{n}} - {f_{wn}\lambda_{n}}}},{B_{an} = {{\omega^{2}\rho_{a}S_{n}\lambda_{n}} - {f_{an}\lambda_{n}}}},{B_{n} = {{2\;\mu\; S_{n}\lambda_{n}^{2}} + \lambda - {a\;\mathcal{X}\; f_{wn}} - {{a\left( {1 - \mathcal{X}} \right)}f_{an}}}},{S_{n} = {- \frac{\lambda + \mu + {b_{1}f_{wn}} + {b_{2}f_{an}}}{b_{3} + {\mu\; d_{n}}}}},{f_{an} = \frac{{b_{13}b_{31}} - {b_{11}b_{33}} + {b_{11}d_{n}}}{{b_{11}b_{32}} - {b_{31}b_{12}} + {b_{31}d_{n}}}},{f_{wn} = \frac{{b_{32}b_{23}} - {b_{22}b_{33}} + {b_{22}d_{n}}}{{b_{22}b_{31}} - {b_{32}b_{21}} + {b_{32}d_{n}}}},{\lambda_{n} = \sqrt{\beta^{2} + \gamma^{2} + d_{n}}},d_{n}$ is roots of an equation

${{d_{n}^{3} + {b_{4}d_{n}^{2}} + {b_{5}d_{n}} + b_{6}} = 0},{b_{4} = {{- b_{12}} - b_{33} - b_{21}}},{b_{5} = {{b_{12}b_{33}} - {b_{13}b_{32}} + {b_{21}b_{12}} + {b_{21}b_{33}} - {b_{22}b_{11}} - {b_{23}b_{31}}}},{b_{6} = {{{- b_{21}}b_{12}b_{33}} - {b_{13}b_{31}b_{22}} - {b_{23}b_{11}b_{32}} + {b_{21}b_{13}b_{32}} + {b_{22}b_{11}b_{33}} + {b_{23}b_{31}b_{12}}}},{n = 1},2,3,\;{b_{1} = {{{- a}\;\mathcal{X}} - \frac{\omega^{2}}{b^{w}}}},{b_{2} = {{- {a\left( {1 - \mathcal{X}} \right)}} - \frac{\omega^{2}}{b^{a}}}},{b_{3} = {{\rho\;\omega^{2}} + \frac{\omega^{4}\rho_{w}}{b^{w}} + \frac{\omega^{4}\rho_{a}}{b_{a}}}},{b_{11} = {b^{a}\rho_{a}A_{11}}},{b_{12} = {b^{a}\rho_{a}A_{12}}},{b_{13} = {{\omega^{2}\rho_{a}} + {b^{a}\rho_{a}A_{13}}}},{b_{21} = {b^{w}\rho_{w}A_{21}}},{b_{22} = {b^{w}\rho_{w}A_{22}}},{b_{23} = {{\omega^{2}\rho_{w}} + {b^{w}\rho_{w}A_{23}}}},{b_{31} = {- \frac{{b_{1}b_{21}} + {b_{2}b_{11}}}{\lambda + {2\;\mu}}}},{b_{32} = {- \frac{{b_{1}b_{22}} + {b_{2}b_{12}}}{\lambda + {2\;\mu}}}},{b_{33} = {- \frac{{b_{1}b_{23}} + {b_{2}b_{13}} + b_{3}}{\lambda + {2\;\mu}}}},{b^{w} = {{- \frac{\omega^{2}}{n\; S_{r}}} + \frac{i\;\omega\; g}{k_{w}}}},{b^{a} = {{- \frac{\omega^{2}}{n\left( {1 - S_{r}} \right)}} + \frac{i\;\omega\; g}{k_{a}}}},$

a total stiffness matrix of a flexible pavement expressway model is: {Σ _(f) }=[K ^(f) ]{Ū _(f)},

wherein {Σ _(f)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ 0 0 0 0 0 0 0 0}^(T), {Ū_(f)}={ū_(x) ⁰ ū_(y) ⁰ iū_(z) ⁰ ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), superscripts 0, 1, 2 correspond to different layer positions shown from top to bottom in the three-layered model, and [K^(f)] is a matrix of 11×11,

a total stiffness matrix of a rigid pavement expressway model is: {Σ _(r) }=[K ^(r) ]{Ū _(r)},

wherein {Σ _(r)}={−τ _(xz,b) ¹ −τ _(yz,b) ¹ −iτ _(zz,b) ¹ 0 0 0 0 0}^(T), {Ū_(r)}={ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), and [K^(r)] is a matrix of 8×8.

The solving unit 304 is configured to determine a load form and a load size according to a specific vehicle, substituting the load form and the load size as boundary conditions into the stiffness matrix of a corresponding model to solve and obtain a solution of a parameter to be solved in a frequency domain, wherein the load form and the load size determined according to the specific vehicle are

${f\left( {x_{1},y_{1},t} \right)} = \left\{ {\begin{matrix} {\frac{{Pe}^{i\Omega_{0}t}}{40I_{1}I_{2}}:\left\{ \begin{matrix} {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + l_{2}} \right)}} \end{matrix} \right.} \\ {0;{others}} \end{matrix},} \right.$

wherein P represents a total axle load of the specific vehicle, and in a moving coordinate system, x=x₁−Vt, y=y₁, z=z₁;

wherein in the transform domain, the load form and the load size determined according to the specific vehicle are

${{\overset{¯}{f}\left( {\beta,\gamma} \right)} = \frac{P\left( {e^{{- i}\gamma d_{w}} + e^{i\gamma d_{w}}} \right)}{10l_{1}l_{2}\beta\gamma}}\text{ }{\left\lbrack {{{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}e^{{- i}\beta d_{c}}} + {{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}\left( {e^{{- i}\beta d_{a}} + e^{{- i}\beta d_{a}}} \right)\left( {e^{{- i}\gamma d_{i}} + e^{i\gamma d_{i}}} \right)}} \right\rbrack,}$

the boundary condition of the rigid pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},b}^{1}\left( {x,y,t} \right)} = {- {\varphi\left( {x,y,t} \right)}}} \\ {{\tau_{{xz},b}^{1}\left( {x,y,t} \right)} = {{\tau_{{yz},b}^{1}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,\ y,t} \right)} = 0}}}} \\ {{u_{{zz},b}^{1}\left( {x,y,t} \right)} = {W\left( {x,y,t} \right)}} \end{matrix},} \right.$

the boundary condition of the rigid pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {- {\overset{¯}{\varphi}\left( {\beta,\gamma} \right)}}} \\ {{{\overset{¯}{\tau}}_{{xz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{{yz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \\ {{{\overset{¯}{u}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {\overset{¯}{W}\left( {\beta,\gamma} \right)}} \end{matrix},} \right.$

wherein

${{\overset{¯}{\varphi}\left( {\beta,\gamma} \right)} = \frac{\overset{¯}{f}\left( {\beta,\gamma} \right)}{1 + {B_{33}{D_{P}\left( {\beta^{2} + \gamma^{2}} \right)}^{2}} - {B_{33}\omega^{2}m_{b}}}},$ and B₃₃ is an element in a third row and a third column of [K^(r)]⁻¹ and corresponds to a vertical displacement of the base layer,

the boundary condition of the flexible pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},p}^{0}\left( {x,y,t} \right)} = {f\left( {x,y,t} \right)}} \\ {{\tau_{xz}^{0}\left( {x,y,t} \right)} = {{\tau_{yz}^{0}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \end{matrix},} \right.$

the boundary condition of the flexible pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},p}^{0}\left( {\beta,\gamma} \right)} = {\overset{¯}{f}\left( {\beta,\gamma} \right)}} \\ {{{\overset{¯}{\tau}}_{xz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{x}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \end{matrix},} \right.$

wherein for the rigid pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the rigid pavement expressway model, and for the flexible pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the flexible pavement expressway model, in such a way to obtain the displacement of each surface, the displacement {ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²} of the top surface of the subgrade is substituted into the stress-strain relationship of the unsaturated subgrade in the transform domain to obtain a stress {τ _(xz,s) τ _(yz,s) τ _(zz,s) p _(w) p _(α)} and a displacement {ū_(x) ² ū_(y) ² ū_(z) ² W ² V ²} at any position inside the subgrade in the frequency domain, and a vibration velocity of v _(ij)=iωū_(ij) and an acceleration of ā_(ij)=−ω²ū_(ij) of the subgrade can be calculated according to the displacement.

The converting unit 305 is configured to use an inverse fast Fourier transform to convert the solution f of the parameter to be solved in the frequency domain into a solution f in the time domain, wherein the solution f in the time domain is a stress {τ_(xz,s) τ_(yz,s) τ_(zz,s) p_(w) p_(α)}, a displacement {u² _(x) u² _(y) u² _(z) W² V²}, a velocity v_(ij) and an acceleration a_(ij).

In an exemplary implementation mode, the apparatus further includes a task acquiring unit 306, which is configured to:

acquire a task of determining the performance of the subgrade of the target expressway section.

In another exemplary implementation mode, the task acquiring unit 306 is specifically configured to:

monitor whether project data on a system are updated; and

if yes, acquire the task of determining the performance of the subgrade of the target expressway section based on the updated target project data.

In another exemplary implementation mode, the task acquiring unit 306 is specifically configured to:

receive the task of determining the performance of the subgrade of the target expressway section as inputted by a user.

In another exemplary implementation mode, the apparatus further includes an outputting unit 307, which is configured to:

output the solution f in the time domain.

In another exemplary implementation mode, the outputting unit 307 is specifically configured to:

output a performance determination report of the subgrade of the target expressway section, wherein the performance determination report includes the solution f in the time domain.

In another exemplary implementation mode, the apparatus further includes a warning unit, which is configured to:

determine whether a value of the solution f in the time domain reaches a predetermined alarm value; and

if yes, output an alarm prompt.

The present disclosure also provides a processing device from the perspective of hardware structure. With reference to FIG. 4 , FIG. 4 shows a structural schematic view of the processing device of the present disclosure. Specifically, the processing device of the present disclosure can comprise a processor 401, a memory 402 and an input/output device 403, the processor 401 is used to realize the step of the method for determining a performance of an unsaturated subgrade of a multi-layered expressway section in e.g. the embodiment corresponding to FIG. 1 when the computer program stored in the memory 402 is executed; or the processor 401 is used to realize the function of each unit in e.g. the embodiment corresponding to FIG. 3 when the computer program stored in the memory 402 is executed, and the memory 402 is used to store the computer program required by the processor 401 for executing the method for determining the performance of the subgrade of the expressway section in the embodiment corresponding to FIG. 1 .

Exemplarily, the computer program can be divided into one or more modules/units that are stored in the memory 402 and executed by the processor 401 to accomplish the present disclosure. The one or more modules/units can be a series of computer program instruction segments capable of achieving particular functions, and the instruction segments are used to describe the execution process of the computer program in a computer device.

The processing device may include, but is not limited to, the processor 401, the memory 402, and the input/output device 403. Those skilled in the art can understand that only the example of the processing device is schematically shown and does not constitute any limitation to the processing device, and the processing device may include more or fewer components than those schematically shown, or combine some components or different components. For example, the processing device can also include a network access device, a bus and the like, and the processor 401, the memory 402 and the input/output device 403, etc. are connected by the bus.

The processor 401 may be a central processing unit (CPU), or other general-purpose processor, digital signal processor (DSP), application specific integrated circuit (ASIC), field-programmable gate array (FPGA) or other programmable logic device, discrete gate, or transistor logic device, discrete hardware component, etc. The general-purpose processor may be a microprocessor or any conventional processor. The processor is the control center of the processing device, which is connected to each part of the entire apparatus by means of various interfaces and lines.

The memory 402 may be used to store a computer program and/or module, and the processor 401 realizes the various functions of the computer device by operating or executing a program and/or module stored in the memory 402, and calling data stored in the memory 402. The memory 402 can mainly include a storage program area and a storage data area, wherein the storage program area may store an operating system, an application program required for realizing at least one function, and the like; the storage data area may store data created according to the use of the processing device. In addition, the memory can include a high-speed random access memory, or non-volatile memory, such as a hard disk, internal memory, plug-in hard disk, smart media card (SMC), secure digital (SD) card, flash card, at least one disk storage device, flash device, or other volatile solid-state storage devices.

When executing a computer program stored in the memory 402, the processor 401 is configured to:

use an Odemark equivalent hypothesis of modulus-thickness to simplify a road structure of a target expressway section into a three-layered model, wherein the three-layered model is

${H_{e} = {\left( \frac{E}{E_{1}} \right)^{1/3}H}},$

wherein H_(e) is an equivalent thickness, E is a modulus of a layer to be transformed, H is a thickness of the layer to be transformed, and E₁ is a modulus of a target layer;

based on the three-layered model, set a base layer of the target expressway section as an elastic medium and set the subgrade of the target expressway section as an unsaturated porous medium to derive and obtain a rigid pavement governing equation, a flexible pavement and base layer governing equation, and an unsaturated subgrade governing equation, wherein the rigid pavement governing equation is:

${{{D_{p}\left\lbrack {\frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial x^{4}} + {2\frac{\partial^{4}{W\left( {x,y,t} \right)}}{{\partial x^{2}}{\partial y^{2}}}} + \frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial y^{4}}} \right\rbrack} + {\rho_{p}h_{p}{\overset{¨}{W}\left( {x,y,t} \right)}} + {\varphi\left( {x,y,t} \right)}} = {f\left( {x,y,t} \right)}},$

wherein D_(p)=E_(p)h_(p) ³/[12(1−μ_(p) ²)] is a bending rigidity of a pavement slab, W(x, y, t) is a vertical displacement of the pavement slab, E_(p), ρ_(p) and h_(p) are respectively a modulus, a density and a thickness of the pavement slab, f(x, y, t) is a traffic load on the pavement slab, and φ(x, y, t) is a reaction force exerted by the base layer on the pavement slab,

the flexible pavement and base layer governing equation is

$\left\{ {\begin{matrix} {{{G_{b}u_{t,{jj}}^{b}} + {\left( {\lambda_{o} + G_{b}} \right)u_{j,{ji}}^{b}}} = {\rho_{b}{\overset{¨}{u}}_{i}^{b}}} \\ {\tau_{ij}^{b} = {{\lambda_{b}\theta_{b}\delta_{ij}} + {2G_{b}\varepsilon_{ij}^{b}}}} \end{matrix},} \right.$

wherein u^(b) is a displacement of an elastic medium layer, G_(b) and λ_(b) are Lame constants of the elastic medium layer, ρ_(b) is a density of the elastic medium layer, τ^(b) is a stress of the elastic medium layer, τ^(b) and ε^(b) respectively represent the stress of the elastic medium layer and a strain of the elastic medium layer, and δ_(ij) is a Kronecker symbol,

the unsaturated subgrade governing equation is

$\left\{ {\begin{matrix} {{- p_{w,i}} = {{\rho_{w}{\overset{¨}{u}}_{i}} + {\frac{\rho_{w}}{nS_{r}}{\overset{¨}{W}}_{i}} + {\frac{\rho_{w}g}{k_{w}}{\overset{.}{W}}_{i}}}} \\ {{- p_{a,i}} = {{\rho_{a}{\overset{¨}{u}}_{i}} + {\frac{\rho_{a}}{n\left( {1 - S_{r}} \right)}{\overset{¨}{V}}_{i}} + {\frac{\rho_{a}g}{k_{a}}{\overset{.}{V}}_{i}}}} \\ {{{A_{11}{\overset{.}{p}}_{w}} + {A_{12}{\overset{.}{p}}_{a}} + {A_{13}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{14}{\nabla \cdot {\overset{.}{V}}_{i}}}} = 0} \\ {{{A_{21}{\overset{.}{p}}_{w}} + {A_{22}{\overset{.}{p}}_{a}} + {A_{23}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{24}{\nabla \cdot {\overset{.}{W}}_{i}}}} = 0} \\ {{{\mu{\nabla^{2}{\overset{.}{u}}_{i}}} + {\left( {\lambda + G} \right){\nabla e}} - {a\chi{\nabla p_{w}}} - {{a\left( {1 - \chi} \right)}{\nabla p_{a}}}} = {{\rho{\overset{¨}{u}}_{i}} + {\rho_{w}{\overset{¨}{W}}_{i}} + {\rho_{a}{\overset{¨}{V}}_{i}}}} \end{matrix},} \right.$

wherein

${A_{11} = {{nA_{s}} + {\left( {1 - S_{r}} \right)\frac{{a\chi} - {nS_{r}}}{K_{s}}}}},$ ${A_{12} = {{\left( {1 - S_{r}} \right)\left\lbrack {\frac{{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}}{K_{s}} + \frac{n}{K_{a}}} \right\rbrack} - {nA_{s}}}},$ ${A_{13} = {\left( {1 - S_{r}} \right)\left( {1 - \frac{K_{b}}{K_{s}}} \right)}},{A_{14} = 1},$ ${A_{21} = {\frac{\left( {{a\chi} - {nS_{r}}} \right)S_{r}}{K_{s}} - {nA_{s}} + \frac{nS_{r}}{K_{w}}}},$ ${A_{22} = {\frac{\left\lbrack {{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}} \right\rbrack S_{r}}{K_{s}} + {nA_{s}}}},$ ${A_{23} = {\left( {1 - \frac{K_{b}}{K_{s}}} \right)S_{r}}},{A_{24} = 1},$ ${A_{s} = {{{- \alpha}{md}\left( 1 \right.} - {\left. S_{w0} \right){\left( S_{e} \right)^{\frac{m + 1}{m}}\left\lbrack {\left( S_{e} \right)^{- \frac{1}{m}} - 1} \right\rbrack}^{\frac{d - 1}{d}}}}},$ a = 1 − K_(b)/K_(s), K_(b)andK_(s) are bulk compressibility moduli of soil skeletons and soil particles respectively, p_(w) and p_(a) are pressures of pore water and gas respectively, ρ_(w) and ρ_(a) are densities of pore water and gas respectively, W_(i) and V_(i) are displacement components of pore water and pore gas relative to soil particles in an i direction, g represents a gravitational acceleration, k_(w) and k_(a) represent permeability coefficients of pore water and pore gas respectively, χ is an effective stress parameter, α, m and d are fitting parameters of a soil-water characteristic curve model, S_(e) is an effective saturation, S_(w0) is an irreducible saturation, n is a soil porosity, and S_(r) is a saturation;

based on an action form of a vehicle load, respectively set time terms of variables in the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation respectively as a simple harmonic form e^(iΩ) ⁰ ^(t) and separate the simple harmonic form e^(iΩ) ⁰ ^(t) of the time terms from the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation, and then use a double Fourier transform to respectively reduce partial derivatives in a three-dimensional space domain of the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as ordinary differential equations, to obtain element stiffness matrices of different layer positions after solving and then establish a total dynamic stiffness matrix, wherein the double Fourier transform is

${{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)} = {\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{f\left( {x,y,z,t} \right)}e^{- {i({{\beta x} + {\gamma y}})}}{dxdy}}}}},$ ${{f\left( {x,y,z,t} \right)} = {\frac{1}{4\pi^{2}}{\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)}e^{i({{\beta x} + {\gamma y}})}d\beta d\gamma}}}}},$

wherein β and γ are wave numbers in x and y directions respectively, and f and f are variables in the space domain and a corresponding transform domain respectively,

the rigid pavement governing equation in the transform domain is: W (β,γ)[D _(p)(β²+γ²)²−Ω²ρ_(p) h _(p)]+φ(β,γ)= f (β,γ),

a stress-strain relationship of the flexible pavement in the transform domain is:

${{{{\left\{ {\sum\limits^{\_}}_{P} \right\}\left\lbrack S^{p} \right\rbrack}\left\lbrack D^{p} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{P} \right\}} = {\left\lbrack K^{p} \right\rbrack\left\{ {\overset{¯}{U}}_{P} \right\}}},$ ${\left\lbrack D^{P} \right\rbrack = \begin{bmatrix} \frac{{- i}\beta{e^{- \alpha_{1}}}^{h}}{k_{1}^{2}} & \frac{{- i}\beta}{k_{1}^{2}} & e^{{- \alpha_{2}}h} & 1 & 0 & 0 \\ \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\gamma}{k_{1}^{2}} & 0 & 0 & e^{{- \alpha_{2}}h} & 1 \\ \frac{{- i}\alpha_{1}{e^{- \alpha_{1}}}^{h}}{k_{1}^{2}} & \frac{i\alpha_{1}}{k_{1}^{2}} & {\frac{\beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \beta}{\alpha_{2}} & {\frac{\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \gamma}{\alpha_{2}} \\ \frac{{- i}\beta}{k_{1}^{2}} & \frac{{- i}\beta e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 1 & e^{{- \alpha_{2}}h} & 0 & 0 \\ \frac{{- i}\gamma}{k_{1}^{2}} & \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 0 & 0 & 1 & e^{{- \alpha_{2}}h} \\ \frac{{- i}\alpha_{1}}{k_{1}^{2}} & \frac{i\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{\beta}{\alpha_{2}} & {\frac{- \beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\gamma}{\alpha_{2}} & {\frac{- \gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \end{bmatrix}},$ ${\left\lbrack S^{P} \right\rbrack = {G_{p}\begin{bmatrix} {\frac{2i{\beta\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i{\beta\alpha}_{1}}{k_{1}^{2}}} & {{- \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha^{2} + \beta^{2}} \right)}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} \\ {\frac{2i\gamma\alpha}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{{- 2}i{\gamma\alpha}_{1}}{k_{1}^{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} \\ {\frac{i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}} & {{- 2}\beta e^{{- \alpha_{2}}h}} & {{- 2}\beta} & {{- 2}\gamma e^{{- \alpha_{2}}h}} & {{- 2}\gamma} \\ \frac{{- 2}i{\beta\alpha}_{1}}{k_{1}^{2}} & {\frac{2i{\beta a}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} \\ \frac{{- 2}i\gamma\alpha_{1}}{k_{1}^{2}} & {\frac{2i{\gamma\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}} & {\frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {2\beta} & {2\beta e^{{- \alpha_{2}}h}} & {2\gamma} & {2\gamma e^{{- \alpha_{2}}h}} \end{bmatrix}}},$

wherein superscript 0 of elements in the matrix represents a top surface of a flexible pavement layer, and superscript 1 of elements in the matrix represents a bottom surface of the flexible pavement layer, a stress vector {Σ _(P)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ τ _(xz,p) ¹ τ _(yz,p) ¹ iτ _(zz,p) ¹}^(T), [D^(p)] is a displacement matrix of the flexible pavement, [S^(p)] is a stress matrix of the flexible pavement, a displacement vector {Ū_(P)}={ū_(x,p) ⁰ ū_(y,p) ⁰ iū_(z,p) ⁰ ū_(x,p) ¹ ū_(y,p) ¹ iū_(z,p) ¹}^(T), [K^(p)]=[S^(p)][D^(p)]⁻¹ is a dynamic stiffness matrix of the flexible pavement layer,

${\alpha_{1}^{2} = {\beta^{2} + \gamma^{2} - k_{1}^{2}}},{\alpha_{2}^{2} = {\beta^{2} + \gamma^{2} - k_{2}^{2}}},^{k_{1}^{2} = \frac{\omega^{2}}{c_{1}^{2}}},^{k_{2}^{2} = \frac{\omega^{2}}{c_{2}^{2}}},$ c is a traffic load velocity, c₁ and c₂ are a compression wave velocity and a shear wave velocity of the elastic medium respectively, ω=Ω₀−βc, h represents a thickness of the flexible pavement layer, and i is an imaginary number,

a stress-strain relationship of the elastic base layer in the transform domain is: {Σ _(b) }=[S ^(b) ][D ^(b)]⁻¹ {Ū _(b) }=[K ^(b) ]{Ū _(b)},

a stress-strain relationship of the unsaturated subgrade in the transform domain is:

${\left\{ {\sum\limits^{\_}}_{s} \right\} = {{{\left\lbrack S^{s} \right\rbrack\left\lbrack D^{s} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{s} \right\}} = {\left\lbrack K^{s} \right\rbrack\left\{ {\overset{¯}{U}}_{s} \right\}}}},$ ${\left\lbrack D^{s} \right\rbrack = \left\lbrack \text{⁠}\begin{matrix} {{- \frac{\gamma}{\beta}}e^{{- \lambda_{0}}z}} & {\frac{i\lambda_{0}}{\beta}e^{{- \lambda_{0}}z}} & {S_{1}i\beta e^{{- \lambda_{1}}z}} & {S_{2}i\beta e^{{- \lambda_{2}}z}} & {S_{3}i\beta e^{{- \lambda_{3}}z}} \\ e^{{- \lambda_{0}}z} & 0 & {S_{1}i\gamma e^{{- \lambda_{1}}z}} & {S_{2}i\gamma e^{{- \lambda_{2}}z}} & {S_{3}i\gamma e^{{- \lambda_{3}}z}} \\ 0 & {{- i}e^{{- \lambda_{0}}z}} & {{- i}S_{1}\lambda_{1}e^{{- \lambda_{1}}z}} & {{- i}S_{2}\lambda_{2}e^{{- \lambda_{2}}z}} & {{- i}S_{3}\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{w}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{w1}}{b^{w}\rho_{w}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{w2}}{b^{w}\rho_{w}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{w3}}{b^{w}\rho_{w}}}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{a}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{a1}}{b^{a}\rho_{a}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{a2}}{b^{a}\rho_{a}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{a3}}{b^{a}\rho_{a}}}e^{{- \lambda_{3}}z}} \end{matrix} \right\rbrack},$

$\left\lbrack S^{s} \right\rbrack = \begin{bmatrix} {{- \frac{{\gamma\lambda}_{0}G}{\beta}}e^{{- \lambda_{0}}z}} & {i{G\left( {\frac{\lambda_{0}^{2}}{\beta} + \beta} \right)}e^{{- \lambda_{0}}z}} & {2GS_{1}i\beta\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i\beta\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i\beta\lambda_{3}e^{{- \lambda_{3}}z}} \\ {G\lambda_{0}e^{{- \lambda_{0}}z}} & {{iG}\gamma e^{{- \lambda_{0}}z}} & {2GS_{1}i{\gamma\lambda}_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i{\gamma\lambda}_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i{\gamma\lambda}_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- i}2G\lambda_{0}e^{{- \lambda_{0}}z}} & {{- i}B_{1}e^{{- \lambda_{1}}z}} & {{- i}B_{2}e^{{- \lambda_{2}}z}} & {{- i}B_{3}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{w1}}e^{{- \lambda_{1}}z}} & {{- f_{w2}}e^{{- \lambda_{2}}z}} & {{- f_{w3}}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{a1}}e^{{- \lambda_{1}}z}} & {{- f_{a2}}e^{{- \lambda_{2}}z}} & {{- f_{a3}}e^{{- \lambda_{3}}z}} \end{bmatrix}$

wherein

$\begin{matrix} {{\left\{ {\overset{\_}{\Sigma}}_{S} \right\} = \left\{ {{- {\overset{¯}{\tau}}_{{xz},s}}\  - {\overset{¯}{\tau}}_{{yz},s}\  - {i{\overset{¯}{\tau}}_{{zz},s}}\  - {\overset{¯}{p}}_{w}\  - {\overset{¯}{p}}_{a}} \right\}^{T}},} \\ {\begin{matrix} {{\left\{ {\overset{¯}{U}}_{s} \right\} = \left\{ {{\overset{¯}{u}}_{x,s}\ {\overset{¯}{u}}_{y,s}\ i{\overset{¯}{u}}_{z,s}\ \overset{¯}{W}\ \overset{¯}{V}} \right\}^{T}},} & {B_{wn} = {{\omega^{2}\rho_{w}S_{n}\lambda_{n}} - {f_{wn}\lambda_{n}}}} \end{matrix},} \\ {\begin{matrix} {{B_{an} = {{\omega^{2}\rho_{w}S_{n}\lambda_{n}} - {f_{an}\lambda_{n}}}},} & {B_{n} = {{2\mu S_{n}\lambda_{n}^{2}} + \lambda - {a\chi f_{wn}} - {a\left( {1 - \chi} \right)f_{an}}}} \end{matrix},} \\ \begin{matrix} {{S_{n} = {- \frac{\lambda + \mu + {b_{1}f_{wn}} + {b_{2}f_{an}}}{b_{3} + {\mu d_{n}}}}},} & {{f_{an} = \frac{{b_{13}b_{31}} - {b_{11}b_{33}} + {b_{11}d_{n}}}{{b_{11}b_{32}} - {b_{31}b_{12}} + {b_{31}d_{n}}}},} & {{f_{wn} = \frac{{b_{32}b_{23}} - {b_{22}b_{33}} + {b_{22}d_{n}}}{{b_{22}b_{31}} - {b_{32}b_{21}} + {b_{32}d_{n}}}},} \end{matrix} \\ {{\lambda_{n} = \sqrt{\beta^{2} + \gamma^{2} + d_{n}}},d_{n}} \end{matrix}$ is roots of an equation

$\begin{matrix} {{{d_{n}^{3} + {b_{4}d_{n}^{2}} + {b_{5}d_{n}} + b_{6}} = 0},{b_{4} = {{- b_{12}} - b_{33} - b_{21}}},} \\ {b_{5} = {{b_{12}b_{33}} - {b_{13}b_{32}} + {b_{21}b_{12}} + {b_{21}b_{33}} - {b_{22}b_{11}} - {b_{23}b_{31}}}} \\ {{b_{6} = {{{- b_{21}}b_{12}b_{33}} - {b_{13}b_{31}b_{22}} - {b_{23}b_{11}b_{32}} + {b_{21}b_{13}b_{32}} + {b_{22}b_{11}b_{33}} + {b_{23}b_{31}b_{12}}}},{n = 1},2,3,{b_{1} = {{- {\alpha\chi}}\frac{\omega^{2}}{b^{w}}}},} \\ {{b_{2} = {{- {a\left( {1 - \chi} \right)}} - \frac{\omega^{2}}{b^{a}}}},{b_{3} = {{\rho\omega}^{2} + \frac{\omega^{4}\rho_{w}}{b^{w}} + \frac{\omega^{4}\rho_{a}}{b^{a}}}},{b_{11} = {b^{a}\rho_{a}A_{11}}},{b_{12} = {b^{a}\rho_{a}A_{12}}},} \\ \begin{matrix} {{b_{13} = {{\omega^{2}\rho_{a}} + {b^{a}\rho_{a}A_{13}}}},} & {{b_{21} = {b^{w}\rho_{w}A_{21}}},} & {{b_{22} = {b^{w}\rho_{w}A_{22}}},} & {{b_{23} = {{\omega^{2}\rho_{w}} + {b^{w}\rho_{w}A_{23}}}},} \end{matrix} \\ {{b_{31} = \frac{{b_{1}b_{21}} + {b_{2}b_{11}}}{\lambda + {2\mu}}},{b_{32} = {- \frac{{b_{1}b_{22}} + {b_{2}b_{12}}}{\lambda + {2\mu}}}},{b_{33} = {- \frac{{b_{1}b_{23}} + {b_{2}b_{13}} + b_{3}}{\lambda + {2\mu}}}},{b^{w} = {{- \frac{\omega^{2}}{nS_{r}}} + \frac{i{\omega g}}{k_{w}}}},} \\ {{b^{a} = {{- \frac{\omega^{2}}{n\left( {1 - S_{r}} \right)}} + \frac{i\omega g}{k_{a}}}},} \end{matrix}$

a total stiffness matrix of a flexible pavement expressway model is: {Σ _(f) }=[K ^(f) ]{Ū _(f)},

wherein {Σ _(f)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ 0 0 0 0 0 0 0 0}^(T), {Ū_(f)}={ū_(x) ⁰ ū_(y) ⁰ iū_(z) ⁰ ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), superscripts 0, 1, 2 correspond to different layer positions shown from top to bottom in the three-layered model, and [K^(f)] is a matrix of 11×11,

a total stiffness matrix of a rigid pavement expressway model is: {Σ _(r) }=[K ^(r) ]{Ū _(r)},

wherein {Σ _(r)}={−τ _(xz,b) ¹ −τ _(yz,b) ¹ −iτ _(zz,b) ¹ 0 0 0 0 0}^(T), {Ū_(r)}={ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), and [K^(r)] is a matrix of 8×8.

The processor 401 is further configured to determine a load form and a load size according to a specific vehicle, substitute the load form and the load size as boundary conditions into the stiffness matrix of a corresponding model to solve and obtain a solution f of the parameter to be solved in a frequency domain, wherein the load form and the load size determined according to the specific vehicle are

${f\left( {x_{1},y_{1},t} \right)} = \left\{ \begin{matrix} {\frac{{Pe}^{i\Omega_{0}t}}{40l_{1}l_{22}};\left\{ \begin{matrix} {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + l_{2}} \right)}} \end{matrix} \right.} \\ {0;{others}} \end{matrix} \right.$

wherein P represents a total axle load of the specific vehicle, and in a moving coordinate system, x=x₁−Vt, y=y₁, z=z₁,

wherein in the transform domain, the load form and the load size determined according to the specific vehicle are

${{\overset{¯}{f}\left( {\beta,\gamma} \right)} = \frac{P\left( {e^{{- i}\gamma d_{w}} + e^{i\gamma d_{w}}} \right)}{10l_{1}l_{2}{\beta\gamma}}}\text{ }{\left\lbrack {{{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}e^{{- i}\beta d_{c}}} + {{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}\left( {e^{{- i}\beta d_{a}} + e^{i\beta d_{a}}} \right)\left( {e^{{- i}\gamma d_{t}} + e^{i\gamma d_{t}}} \right)}} \right\rbrack,}$

a boundary condition of the rigid pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},b}^{1}\left( {x,y,t} \right)} = {- {\varphi\left( {x,y,t} \right)}}} \\ {{\tau_{{xz},b}^{1}\left( {x,y,t} \right)} = {{\tau_{{yz},b}^{1}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \\ {{u_{{zz},b}^{1}\left( {x,y,t} \right)} = {W\left( {x,y,t} \right)}} \end{matrix},} \right.$

the boundary condition of the rigid pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {- {\overset{¯}{\varphi}\left( {\beta,\gamma} \right)}}} \\ {{{\overset{¯}{\tau}}_{{xz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{{yz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \\ {{{\overset{¯}{u}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {\overset{¯}{W}\left( {\beta,\gamma} \right)}} \end{matrix},} \right.$

wherein

${{\overset{\_}{\varphi}\left( {\beta,\gamma} \right)} = \frac{\overset{¯}{f}\left( {\beta,\gamma} \right)}{1 + {B_{33}{D_{P}\left( {\beta^{2} + \gamma^{2}} \right)}^{2}} - {B_{33}\omega^{2}m_{b}}}},$ and B₃₃ is an element in a third row and a third column of [K^(r)]⁻¹ and corresponds to a vertical displacement of the base layer,

a boundary condition of the flexible pavement expressway model is:

$\left\{ {\begin{matrix} {{\tau_{{zz},p}^{0}\left( {x,y,t} \right)} = {f\left( {x,y,t} \right)}} \\ {{\tau_{xz}^{0}\left( {x,y,t} \right)} = {{\tau_{yz}^{0}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \end{matrix},} \right.$

the boundary condition of the flexible pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as:

$\left\{ {\begin{matrix} {{{\overset{¯}{\tau}}_{{zz},p}^{0}\left( {\beta,\gamma} \right)} = {\overset{¯}{f}\left( {\beta,\gamma} \right)}} \\ {{{\overset{¯}{\tau}}_{xz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{yz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \end{matrix},} \right.$

wherein for the rigid pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the rigid pavement expressway model, and for the flexible pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the flexible pavement expressway model, in such a way to obtain the displacement of each surface, the displacement {ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²} of the top surface of the subgrade is substituted into the stress-strain relationship of the unsaturated subgrade in the transform domain to obtain a stress {τ _(xz,s) τ _(yz,s) τ _(zz,s) p _(w) p _(α)} and a displacement {ū_(x) ² ū_(y) ² ū_(z) ² W ² V ²} at any position inside the subgrade in the frequency domain, and the vibration velocity of v _(ij)=iωū_(ij) and an acceleration of ā_(ij)=ω²ū_(ij) of the subgrade can be calculated according to the displacement.

The processor 401 is further configured to us an inverse fast Fourier transform to convert the solution f of the parameter to be solved in the frequency domain into a solution f in the time domain, wherein the solution f in the time domain is a stress {τ_(xz,s) τ_(yz,s) τ_(zz,s) p_(w) p_(α)}, a displacement {u² _(x) u² _(y) u² _(z) W² V₂}, a velocity v_(ij) and an acceleration a_(ij).

Those skilled in the art can clearly understand that, for the convenience and brevity of the description, regarding the apparatus for determining the performance of the subgrade of the expressway section, a processing device, and the specific working process of a corresponding unit thereof as described above, reference can be made to the explanation of the method for determining the performance of the subgrade of the expressway section in e.g. the embodiment corresponding to FIG. 1 , which will not be reiterated herein.

Those ordinarily skilled in the art can understand that all or part of the steps of the various methods in the above embodiments can be completed by instructions, or by relevant hardware controlled by instructions. The instructions can be stored in a computer readable storage medium, and loaded and executed by a processor.

Therefore, the present disclosure provides a computer readable storage medium, in which a plurality of instructions is stored. The instructions can be loaded by the processor to execute the steps of the method for determining a performance of an unsaturated subgrade of a multi-layered expressway section in e.g. the embodiment corresponding to FIG. 1 of the present disclosure. For specific operation, reference can be made to the explanation of the method for determining the performance of the subgrade of the expressway section in e.g. the embodiment corresponding to FIG. 1 , which will not be reiterated herein.

Wherein the computer readable storage medium may include a read-only memory (ROM), a random access memory (RAM), a disk or a compact disk, etc.

Since the instruction stored in the computer readable storage medium can execute the steps of the method in e.g. the embodiment corresponding to FIG. 1 of the present disclosure, the advantageous effects of the method for determining a performance of an unsaturated subgrade of a multi-layered expressway section in e.g. the embodiment corresponding to FIG. 1 of the present disclosure can be achieved as stated above in detail, which will not be reiterated herein.

The method and apparatus for determining a performance of an unsaturated subgrade of a multi-layered expressway section, the processing device and the computer readable storage medium provided by the present disclosure are introduced above in detail. Specific examples are used herein to expound the principle and implementation modes of the present disclosure. The above explanations of embodiments are only used to help people understand the method and core idea of the present disclosure; and meanwhile, for those skilled in the art and according to the idea of the present disclosure, there may be changes in the specific implementation modes and application scope. In summary, the contents of the specification should not be construed as any limitation to the present disclosure. 

What is claimed is:
 1. A method for determining a performance of a subgrade of an expressway section, comprising: Step S101: a processing device using an Odemark equivalent hypothesis of modulus-thickness to simplify a road structure of a target expressway section into a three-layered model; Step S102: based on the three-layered model, the processing device setting a base layer of the target expressway section as an elastic medium and setting the subgrade of the target expressway section as an unsaturated porous medium, to derive and obtain a rigid pavement governing equation, a flexible pavement and base layer governing equation and an unsaturated subgrade governing equation; Step S103: based on an action form of a vehicle load, the processing device respectively setting time terms of variables in the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as a simple harmonic form e^(iΩ) ⁰ ^(t) and separating the simple harmonic form e^(iΩ) ⁰ ^(t) of the time terms from the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation, and then using a double Fourier transform to respectively reduce partial derivatives in a three-dimensional space domain of the rigid pavement governing equation, the flexible pavement and base layer governing equation and the unsaturated subgrade governing equation as ordinary differential equations, to obtain element stiffness matrices of different layer positions after solving and then establish a total dynamic stiffness matrix; Step S104: the processing device determining a load form and a load size according to a specific vehicle, substituting the load form and the load size as boundary conditions into the total dynamic stiffness matrix of a corresponding model to solve and obtain a solution f of a parameter to be solved in a frequency domain; and Step S105: the processing device using an inverse fast Fourier transform to convert the solution f of the parameter to be solved in the frequency domain into a solution f in a time domain, wherein the solution fin the time domain is a stress {τ_(xz,s) τ_(yz,s) τ_(zz,s) p_(w) p_(α)}, a displacement {u² _(x) u² _(y) u² _(z) W² V²}, a velocity v_(ij) and an acceleration a_(ij); wherein, in the Step S101, the three-layered model is: ${H_{e} = {\left( \frac{E}{E_{1}} \right)^{1/3}H}},$  wherein H_(e) is an equivalent thickness, E is a modulus of a layer to be transformed, H is a thickness of the layer to be transformed, and E₁ is a modulus of a target layer; wherein, in the step S102, the rigid pavement governing equation is: ${{D_{p}\left\lbrack {\frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial x^{4}} + {2\frac{\partial^{4}{W\left( {x,y,t} \right)}}{{\partial x^{2}}{\partial y^{2}}}} + \frac{\partial^{4}{W\left( {x,y,t} \right)}}{\partial y^{4}}} \right\rbrack} + {\rho_{p}h_{p}{\overset{¨}{W}\left( {x,y,t} \right)}} + {\varphi\left( {x,y,t} \right)}} = {f\left( {x,y,t} \right)}$ wherein D_(p)=E_(p)h_(p) ³/[12(1−μ_(p) ²)] is a bending rigidity of a pavement slab, W(x, y, t) is a vertical displacement of the pavement slab, E_(p), ρ_(p) and h_(p) are respectively a modulus, a density and a thickness of the pavement slab, f(x, y, t) is a traffic load on the pavement slab, and φ(x, y, t) is a reaction force exerted by the base layer on the pavement slab, the flexible pavement and base layer governing equation is: $\left\{ {\begin{matrix} {{{G_{b}u_{i,{jj}}^{b}} + {\left( {\lambda_{b} + G_{b}} \right)u_{j,{ji}}^{b}}} = {\rho_{b}{\overset{¨}{u}}_{i}^{b}}} \\ {\tau_{ij}^{b} = {{\lambda_{b}\theta_{b}\delta_{ij}} + {2G_{b}\varepsilon_{ij}^{b}}}} \end{matrix},} \right.$ wherein u^(b) is a displacement of an elastic medium layer, G_(b) and λ_(b) are Lame constants of the elastic medium layer, ρ_(b) is a density of the elastic medium layer, τ^(b) is a stress of the elastic medium layer, τ^(b) and ε^(b) respectively represent the stress of the elastic medium layer and a strain of the elastic medium layer, and δ_(ij) is a Kronecker symbol, wherein, the unsaturated subgrade governing equation is: $\left\{ \begin{matrix} {{- p_{w,t}} = {{\rho_{w}{\overset{¨}{u}}_{i}} + {\frac{\rho_{w}}{nS_{r}}{\overset{¨}{W}}_{i}} + {\frac{\rho_{w}g}{k_{w}}{\overset{.}{W}}_{i}}}} \\ {{- p_{a,i}} = {{\rho_{a}{\overset{¨}{u}}_{i}} + {\frac{\rho_{a}}{n\left( {1 - S_{r}} \right)}{\overset{¨}{V}}_{i}} + {\frac{\rho_{a}g}{k_{a}}{\overset{.}{V}}_{i}}}} \\ {{{A_{11}{\overset{.}{p}}_{w}} + {A_{12}{\overset{.}{p}}_{a}} + {A_{13}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{14}{\nabla \cdot {\overset{.}{V}}_{i}}}} = 0} \\ {{{A_{21}{\overset{.}{p}}_{w}} + {A_{22}{\overset{.}{p}}_{a}} + {A_{23}{\nabla \cdot {\overset{.}{u}}_{i}}} + {A_{24}{\nabla \cdot {\overset{.}{W}}_{i}}}} = 0} \\ {{{\mu{\nabla^{2}{\overset{.}{\mu}}_{i}}} + {\left( {\lambda + G} \right){\nabla e}} - {a\chi{\nabla p_{w}}} - {{a\left( {1 - \chi} \right)}{\nabla p_{a}}}} = {{\rho{\overset{¨}{u}}_{i}} + {\rho_{w}{\overset{¨}{W}}_{i}} + {\rho_{a}{\overset{¨}{V}}_{i}}}} \end{matrix} \right.$ wherein ${A_{11} = {{nA_{s}} + {\left( {1 - S_{r}} \right)\frac{{a\chi} - {nS}_{r}}{K_{s}}}}},$ ${A_{12} = {{\left( {1 - S_{r}} \right)\left\lbrack {\frac{{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}}{K_{s}} + \frac{n}{K_{a}}} \right\rbrack} - {nA_{s}}}},$ ${A_{13} = {\left( {1 - S_{r}} \right)\left( {1 - \frac{K_{b}}{K_{s}}} \right)}},$ A₁₄ = 1, ${A_{21} = {\frac{\left( {{a\chi} - {nS_{r}}} \right)S_{r}}{K_{s}} - {nA_{s}} + \frac{nS_{r}}{K_{w}}}},$ ${A_{22} = {\frac{\left\lbrack {{a\left( {1 - \chi} \right)} - {n\left( {1 - S_{r}} \right)}} \right\rbrack S_{r}}{K_{s}} + {nA_{s}}}},$ ${A_{23} = {\left( {1 - \frac{K_{b}}{K_{s}}} \right)S_{r}}},$ A₂₄ = 1, ${A_{s} = {{- {{\alpha{md}}\left( {1 - S_{w0}} \right)}}{\left( S_{e} \right)^{\frac{m + 1}{m}}\left\lbrack {\left( S_{e} \right)^{- \frac{1}{m}} - 1} \right\rbrack}^{\frac{d - 1}{d}}}},$ a = 1 − K_(b)/K_(s), K_(b) and K_(s) are bulk compressibility moduli of soil skeletons and soil particles respectively, p_(w) and p_(a) are pressures of pore water and gas respectively, ρ_(w) and ρ_(a) are densities of pore water and gas respectively, W_(i) and V_(i) are displacement components of pore water and pore gas relative to soil particles in an i direction, g represents a gravitational acceleration, k_(w) and k_(a) represent permeability coefficients of pore water and pore gas respectively, χ is an effective stress parameter, α, m and d are fitting parameters of a soil-water characteristic curve model, S_(e) is an effective saturation, S_(w0) is an irreducible saturation, n is a soil porosity, and S_(r) is a saturation; wherein, in the Step S103, the double Fourier transform is: ${{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)} = {\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{f\left( {x,y,z,t} \right)}e^{- {i({{\beta x} + {\gamma y}})}}{dxdy}}}}},$ ${{f\left( {x,y,z,t} \right)} = {\frac{1}{4\pi^{2}}{\int_{- \infty}^{\infty}{\int_{- \infty}^{\infty}{{\overset{¯}{f}\left( {\beta,\gamma,z,t} \right)}e^{i({{\beta x} + {\gamma y}})}d\beta d\gamma}}}}},$ wherein β and γ are wave numbers in x and y directions respectively, and f and f are variables in the three-dimensional space domain and a corresponding transform domain respectively, the rigid pavement governing equation in the transform domain is: W (β,γ)[D _(p)(β²+γ²)²−Ω²ρ_(p) h _(p)]+φ(β,γ)= f (β,γ), a stress-strain relationship of the flexible pavement in the transform domain is: ${\left\{ {\sum\limits^{\_}}_{P} \right\} = {{{\left\lbrack S^{p} \right\rbrack\left\lbrack D^{p} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{P} \right\}} = {\left\lbrack K^{p} \right\rbrack\left\{ {\overset{¯}{U}}_{P} \right\}}}},$ ${\left\lbrack D^{P} \right\rbrack = \text{ }\begin{bmatrix} \frac{{- i}\beta e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\beta}{k_{1}^{2}} & e^{{- \alpha_{2}}h} & 1 & 0 & 0 \\ \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{{- i}\gamma}{k_{1}^{2}} & 0 & 0 & e^{{- \alpha_{2}}h} & 1 \\ \frac{{- i}\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{i\alpha_{1}}{k_{1}^{2}} & {\frac{\beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \beta}{\alpha_{2}} & {\frac{\gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{- \gamma}{\alpha_{2}} \\ \frac{{- i}\beta}{k_{1}^{2}} & \frac{{- i}\beta e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 1 & e^{{- \alpha_{2}}h} & 0 & 0 \\ \frac{{- i}\gamma}{k_{1}^{2}} & \frac{{- i}\gamma e^{{- \alpha_{1}}h}}{k_{1}^{2}} & 0 & 0 & 1 & e^{{- \alpha_{2}}h} \\ \frac{{- i}\alpha_{1}}{k_{1}^{2}} & \frac{i\alpha_{1}e^{{- \alpha_{1}}h}}{k_{1}^{2}} & \frac{\beta}{\alpha_{2}} & {\frac{- \beta}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\gamma}{\alpha_{2}} & {\frac{- \gamma}{\alpha_{2}}e^{{- \alpha_{2}}h}} \end{bmatrix}},$ ${\left\lbrack S^{P} \right\rbrack = {G_{p}\begin{bmatrix} {\frac{2i\beta\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i\beta\alpha_{1}}{k_{1}^{2}}} & {{- \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} \\ {\frac{2i\gamma\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {- \frac{2i\gamma\alpha_{1}}{k_{1}^{2}}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} \\ {\frac{i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}{k_{1}^{2}} & {{- 2}\beta e^{{- \alpha_{2}}h}} & {{- 2}\beta} & {{- 2}\gamma e^{{- \alpha_{2}}h}} & {{- 2}\gamma} \\ \frac{{- 2}i\beta\alpha_{1}}{k_{1}^{2}} & {\frac{2i{\beta\alpha}_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \beta^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} \\ \frac{{- 2}i\gamma\alpha_{1}}{k_{1}^{2}} & {\frac{2i\gamma\alpha_{1}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & \frac{\beta\gamma}{\alpha_{2}} & {{- \frac{\beta\gamma}{\alpha_{2}}}e^{{- \alpha_{2}}h}} & \frac{\left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}} & {\frac{- \left( {\alpha_{2}^{2} + \gamma^{2}} \right)}{\alpha_{2}}e^{{- \alpha_{2}}h}} \\ \frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}} & {\frac{- {i\left( {{2\alpha_{2}^{2}} + k_{2}^{2}} \right)}}{k_{1}^{2}}e^{{- \alpha_{1}}h}} & {2\beta} & {2\beta e^{{- \alpha_{2}}h}} & {2\gamma} & {2\gamma e^{{- \alpha_{2}}h}} \end{bmatrix}}},$ wherein a stress vector {Σ _(P)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ τ _(xz,p) ¹ −τ _(yz,p) ¹ iτ _(zz,p) ¹}^(T), [D^(p)] is a displacement matrix of the flexible pavement, [S^(p)] is a stress matrix of the flexible pavement, a displacement vector {Ū_(P)}={ū_(x,p) ⁰ ū_(y,p) ⁰ iū_(z,p) ⁰ ū_(x,p) ¹ ū_(y,p) ¹ iū_(z,p) ¹}^(T), [K^(p)]=[S^(p)][D^(p)]⁻¹ is a dynamic stiffness matrix of the flexible pavement layer, α₁² = β² + γ² − k₁², α₂² = β² + γ² − k₂², ${k_{1}^{2} = \frac{\omega^{2}}{c_{1}^{2}}},$ ${k_{2}^{2} = \frac{\omega^{2}}{c_{2}^{2}}},$  c is a traffic load velocity, c₁ and c₂ are a compression wave velocity and a shear wave velocity of the elastic medium, respectively, ω=Ω₀−βc, h represents a thickness of the flexible pavement layer, and i is an imaginary number, wherein superscript 0 of elements in the stress vector and in the displacement vector represents a top surface of a flexible pavement layer, and superscript 1 of elements in the stress vector and in the displacement vector represents a bottom surface of the flexible pavement layer; wherein a stress-strain relationship of the elastic base layer in the transform domain is: {Σ _(b) }=[S ^(b) ][D ^(b)]⁻¹ {Ū _(b) }=[K ^(b) ]{Ū _(b)}, a stress-strain relationship of the unsaturated subgrade in the transform domain is: ${\left\{ {\sum\limits^{\_}}_{s} \right\} = {{{\left\lbrack S^{s} \right\rbrack\left\lbrack D^{s} \right\rbrack}^{- 1}\left\{ {\overset{¯}{U}}_{s} \right\}} = {\left\lbrack K^{s} \right\rbrack\left\{ {\overset{¯}{U}}_{s} \right\}}}},$ ${\left\lbrack D^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma}{\beta}}e^{- \lambda_{0^{z}}}} & {\frac{i\lambda_{0}}{\beta}e^{- \lambda_{0^{z}}}} & {S_{1}i\beta e^{{- \lambda_{1}}z}} & {S_{2}i\beta e^{- \lambda_{2^{z}}}} & {S_{3}i\beta e^{{- \lambda_{3}}z}} \\ e^{- \lambda_{0^{z}}} & 0 & {S_{1}i\gamma e^{{- \lambda_{1}}z}} & {S_{2}i\gamma e^{- \lambda_{2^{z}}}} & {S_{3}i\gamma e^{{- \lambda_{3}}z}} \\ 0 & {{- i}e^{- \lambda_{0^{z}}}} & {{- i}S_{1}\lambda_{1}e^{{- \lambda_{1}}z}} & {{- i}S_{2}\lambda_{2}e^{{- \lambda_{2}}z}} & {{- i}S_{3}\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{w}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{w1}}{b^{w}\rho_{w}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{w2}}{b^{w}\rho_{w}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{w3}}{b^{w}\rho_{w}}}e^{{- \lambda_{3}}z}} \\ 0 & {{- \frac{\omega^{2}}{b^{a}}}e^{{- \lambda_{0}}z}} & {{- \frac{B_{a1}}{b^{a}\rho_{a}}}e^{{- \lambda_{1}}z}} & {{- \frac{B_{a2}}{b^{a}\rho_{a}}}e^{{- \lambda_{2}}z}} & {{- \frac{B_{a3}}{b^{a}\rho_{a}}}e^{{- \lambda_{3}}z}} \end{bmatrix}_{▯}}$ ${\left\lbrack S^{s} \right\rbrack = \begin{bmatrix} {{- \frac{\gamma\lambda_{0}G}{\beta}}e^{- \lambda_{0^{z}}}} & {{{iG}\left( {\frac{\lambda_{0}^{2}}{\beta} + \beta} \right)}e^{- \lambda_{0^{z}}}} & {2GS_{1}i\beta\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i\beta\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i\beta\lambda_{3}e^{{- \lambda_{3}}z}} \\ {G\lambda_{0}e^{- \lambda_{0^{z}}}} & {{iG}\gamma e^{- \lambda_{0^{z}}}} & {2GS_{1}i\gamma\lambda_{1}e^{{- \lambda_{1}}z}} & {2GS_{2}i\gamma\lambda_{2}e^{{- \lambda_{2}}z}} & {2GS_{3}i\gamma\lambda_{3}e^{{- \lambda_{3}}z}} \\ 0 & {{- i}2G\lambda_{0}e^{- \lambda_{0^{z}}}} & {{- i}B_{1}e^{{- \lambda_{1}}z}} & {{- i}B_{2}e^{{- \lambda_{2}}z}} & {{- i}B_{3}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{w1}}e^{{- \lambda_{1}}z}} & {{- f_{w2}}e^{{- \lambda_{2}}z}} & {{- f_{w3}}e^{{- \lambda_{3}}z}} \\ 0 & 0 & {{- f_{a1}}e^{{- \lambda_{1}}z}} & {{- f_{a2}}e^{{- \lambda_{2}}z}} & {{- f_{a3}}e^{{- \lambda_{3}}z}} \end{bmatrix}},$ wherein $\left\{ {\sum\limits^{\_}}_{s} \right\} = \left\{ {\begin{matrix} {- {\overset{¯}{\tau}}_{{xz},s}} & {- {\overset{¯}{\tau}}_{{yz},s}} & {{- i}{\overset{¯}{\tau}}_{{zz},s}} & {- {\overset{¯}{p}}_{w}} & \left. {- {\overset{¯}{p}}_{a}} \right\}^{T} \end{matrix},} \right.$ $\left\{ {\overset{\_}{U}}_{s} \right\} = \left\{ {\begin{matrix} {\overset{¯}{u}}_{x,s} & {\overset{¯}{u}}_{y,s} & {i{\overset{¯}{u}}_{z,s}} & \overset{¯}{W} & \left. \overset{¯}{V} \right\}^{T} \end{matrix},} \right.$ B_(wn) = ω²ρ_(w)S_(n)λ_(n) − f_(wn)λ_(n), B_(an) = ω²ρ_(a)S_(n)λ_(n) − f_(an)λ_(n), B_(n) = 2μS_(n)λ_(n)² + λ − aχf_(wn) − a(1 − χ)f_(an), ${S_{n} = {- \frac{\lambda + \mu + {b_{1}f_{wn}} + {b_{2}f_{an}}}{b_{3} + {\mu d_{n}}}}},$ ${f_{an} = \frac{{b_{13}b_{31}} - {b_{11}b_{33}} + {b_{11}d_{n}}}{{b_{11}b_{32}} - {b_{31}b_{12}} + {b_{31}d_{n}}}},$ ${f_{wn} = \frac{{b_{32}b_{23}} - {b_{22}b_{33}} + {b_{22}d_{n}}}{{b_{22}b_{31}} - {b_{32}b_{21}} + {b_{32}d_{n}}}},$ ${\lambda_{n} = \sqrt{\beta^{2} + \gamma^{2} + d_{n}}},$  d_(n) is roots of an equation d_(n)³ + b₄d_(n)² + b₅d_(n) + b₆ = 0, b₄ = −b₁₂ − b₃₃ − b₂₁, b₅ = b₁₂b₃₃ − b₁₃b₃₂ + b₂₁b₁₂ + b₂₁b₃₃ − b₂₂b₁₁ − b₂₃b₃₁, b₆ = −b₂₁b₁₂b₃₃ − b₁₃b₃₁b₂₂ − b₂₃b₁₁b₃₂ + b₂₁b₁₃b₃₂ + b₂₂b₁₁b₃₃ + b₂₃b₃₁b₁₂, ${n = 1},2,3,{b_{1} = {{{- a}\chi} - \frac{\omega^{2}}{b^{w}}}},$ ${b_{2} = {{- {a\left( {1 - \chi} \right)}} - \frac{\omega^{2}}{b^{a}}}},$ ${b_{3} = {{\rho\omega}^{2} + \frac{\omega^{4}\rho_{w}}{b^{w}} + \frac{\omega^{4}\rho_{a}}{b^{a}}}},$ b₁₁ = b^(a)ρ_(a)A₁₁, b₁₂ = b^(a)ρ_(a)A₁₂, b₁₃ = ω²ρ_(a) + b^(a)ρ_(a)A₁₃, b₂₁ = b^(w)ρ_(w)A₂₁, b₂₂ = b^(w)ρ_(w)A₂₂, b₂₃ = ω²ρ_(w) + b^(w)ρ_(w)A₂₃, ${b_{31} = {- \frac{{b_{1}b_{21}} + {b_{2}b_{11}}}{\lambda + {2\mu}}}},$ ${{b_{32} = {- \frac{{b_{1}b_{22}} + {b_{2}b_{12}}}{\lambda + {2\mu}}}},}$ ${b_{33} = {- \frac{{b_{1}b_{23}} + {b_{2}b_{13}} + b_{3}}{\lambda + {2\mu}}}},$ ${b^{w} = {{- \frac{\omega^{2}}{nS_{r}}} + \frac{i\omega g}{k_{w}}}},$ ${b^{a} = {{- \frac{\omega^{2}}{n\left( {1 - S_{r}} \right)}} + \frac{i\omega g}{k_{a}}}},$ a total stiffness matrix of a flexible pavement expressway model is: {Σ _(f) }=[K ^(f) ]{Ū _(f)}, wherein {Σ _(f)}={−τ _(xz,p) ⁰ −τ _(yz,p) ⁰ −iτ _(zz,p) ⁰ 0 0 0 0 0 0 0 0}^(T), {Ū_(f)}={ū_(x) ⁰ ū_(y) ⁰ iū_(z) ⁰ ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), superscripts 0, 1, 2 correspond to different layer positions shown from top to bottom in the three-layered model, and [K^(f)] is a matrix of 11×11, a total stiffness matrix of a rigid pavement expressway model is: {Σ _(r) }=[K ^(r) ]{Ū _(r)}, wherein {Σ _(r)}={−τ _(xz,b) ¹ −τ _(yz,b) ¹ −iτ _(zz,b) ¹ 0 0 0 0 0}^(T), {Ū_(r)}={ū_(x) ¹ ū_(y) ¹ iū_(z) ¹ ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²}^(T), and [K^(r)] is a matrix of 8×8; wherein in the Step S104, the load form and the load size determined according to the specific vehicle are: ${f\left( {x_{1},y_{1},t} \right)} = \text{ }\left\{ \begin{matrix} {\frac{{Pe}^{i\Omega_{0}t}}{40l_{1}l_{2}};\left\{ {\begin{matrix} {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {d_{a} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{a} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} + d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + d_{t} + l_{2}} \right)}} \\ {{\left( {{- d_{a}} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {{- d_{a}} + l_{1}} \right)},{\left( {d_{w} - d_{t} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} - d_{t} + l_{2}} \right)}} \\ {{\left( {d_{c} - l_{1}} \right) \leq {x_{1} - {Vt}} \leq \left( {d_{c} + l_{1}} \right)},{\left( {d_{w} - l_{2}} \right) \leq {❘y_{1}❘} \leq \left( {d_{w} + l_{2}} \right)}} \end{matrix},} \right.} \\ {0;{others}} \end{matrix} \right.$ wherein P represents a total axle load of the specific vehicle, and in a moving coordinate system, x=x₁−Vt, y=y₁, z=z₁, wherein in the transform domain, the load form and the load size determined based on the specific vehicle are: ${{\overset{¯}{f}\left( {\beta,\gamma} \right)} = {\frac{P\left( {e^{{- i}\gamma d_{w}} + e^{i\gamma d_{w}}} \right)}{10l_{1}l_{2}\beta\gamma}\left\lbrack \text{⁠}{{{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}e^{{- i}\beta d_{c}}} + {{\sin\left( {\beta l_{1}} \right)}{\sin\left( {\gamma l_{2}} \right)}\left( {e^{{- i}\beta d_{a}} + e^{i\beta d_{a}}} \right)\left( {e^{{- i}\gamma d_{t}} + e^{i\gamma d_{t}}} \right)}} \right\rbrack}},$ a boundary condition of the rigid pavement expressway model is: $\left\{ {\begin{matrix} {{\tau_{{zz},b}^{1}\left( {x,y,t} \right)} = {- {\varphi\left( {x,y,t} \right)}}} \\ {{\tau_{{xz},b}^{1}\left( {x,y,t} \right)} = {{\tau_{{yz},b}^{1}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \\ {{u_{{zz},b}^{1}\left( {x,y,t} \right)} = {W\left( {x,y,t} \right)}} \end{matrix},} \right.$ the boundary condition of the rigid pavement expressway model with a time term e^(iΩ) ⁰ ^(t) removed is expressed as: $\left\{ {\begin{matrix} {{{\overset{\_}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {- {\overset{¯}{\varphi}\left( {\beta,\gamma} \right)}}} \\ {{{\overset{\_}{\tau}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{\tau}}_{{yz},b}^{1}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \\ {{{\overset{\_}{u}}_{{zz},b}^{1}\left( {\beta,\gamma} \right)} = {\overset{\_}{W}\left( {\beta,\gamma} \right)}} \end{matrix},} \right.$ wherein ${{\overset{¯}{\varphi}\left( {\beta,\gamma} \right)} = \frac{\overset{¯}{f}\left( {\beta,\gamma} \right)}{1 + {B_{33}{D_{P}\left( {\beta^{2} + \gamma^{2}} \right)}^{2}} - {B_{33}\omega^{2}m_{b}}}},$  and B₃₃ is an element in a third row and a third column of [K^(r)]⁻¹ and corresponds to a vertical displacement of the base layer, a boundary condition of the flexible pavement expressway model is: $\left\{ {\begin{matrix} {{\tau_{{zz},p}^{0}\left( {x,y,t} \right)} = {f\left( {x,y,t} \right)}} \\ {{\tau_{xz}^{0}\left( {x,y,t} \right)} = {{\tau_{yz}^{0}\left( {x,y,t} \right)} = {{p_{w}^{2}\left( {x,y,t} \right)} = {{p_{a}^{2}\left( {x,y,t} \right)} = 0}}}} \end{matrix},} \right.$ the boundary condition of the flexible pavement expressway model with the time term e^(iΩ) ⁰ ^(t) removed is expressed as: $\left\{ {\begin{matrix} {{{\overset{\_}{\tau}}_{{zz},p}^{0}\left( {\beta,\gamma} \right)} = {\overset{¯}{f}\left( {\beta,\gamma} \right)}} \\ {{{\overset{\_}{\tau}}_{xz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{\tau}}_{yz}^{0}\left( {\beta,\gamma} \right)} = {{{\overset{\_}{p}}_{w}^{2}\left( {\beta,\gamma} \right)} = {{{\overset{¯}{p}}_{a}^{2}\left( {\beta,\gamma} \right)} = 0}}}} \end{matrix},} \right.$ wherein for the rigid pavement expressway model, a corresponding boundary condition is substituted into the total stiffness matrix of the rigid pavement expressway model, and for the flexible pavement expressway model, the corresponding boundary condition is substituted into the total stiffness matrix of the flexible pavement expressway model, in such a way to obtain a displacement of each surface, then a displacement {ū_(x) ² ū_(y) ² iū_(z) ² W ² V ²} of the top surface of the subgrade is substituted into the stress-strain relationship of the unsaturated subgrade in the transform domain to obtain the a stress {τ _(xz,s) τ _(yz,s) τ _(zz,s) p _(w) p _(α)} and a displacement {ū_(x) ² ū_(y) ² ū_(z) ² W ² V ²} at any position inside the subgrade in the frequency domain, and a vibration velocity v _(ij)=iωū_(ij) and an acceleration ā_(ij)=−ω²ū_(ij) of the subgrade can be calculated according to the displacement.
 2. The method according to claim 1, wherein prior to the Step S101, the method further comprises: Step S106: the processing device acquiring a task of determining the performance of the subgrade of the target expressway section.
 3. The method according to claim 1, wherein subsequent to the Step S105, the method further comprises: Step S107: the processing device outputting the solution fin the time domain.
 4. A processing device, comprising a processor and a memory in which a computer program is stored, wherein when the computer program in the memory is executed by the processor, the processor is configured to execute the method according to claim
 1. 5. A non-transitory computer readable storage medium which stores a plurality of instructions that are suitable to be loaded by a processor to execute the method according to claim
 1. 